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Abstract 

The  purpose  of  this  research  was  to  present  the  method  of  Linear  Regression  as  a 
parameter  identification  method  to  determine  the  longitudinal  dimensional  stability 
derivatives  of  a  tactical  missile.  Missile  flight  histories  are  characterized  by  rapid 
accelerations,  rapidly  changing  mass  property  characteristics  with  often  short  flight 
times.  These  characteristics  make  accurate  parameter  estimation  of  the  missile 
aerodynamics  more  challenging  than  for  aircraft.  The  simulation  used  for  this  research 
was  created  in  MATLAB/SIMULINK  based  on  the  missile  trajectory  program,  TRAP. 
The  aerodynamic  data  for  the  6-DoF  missile  model  was  based  on  a  supersonic,  tail 
controlled  missile  similar  to  an  AIM-9X  missile.  Two  command  input  types  were 
investigated  to  determine  if  either  could  induce  an  excitation  of  the  system  modes  of  the 
plant  being  measured  to  lead  to  good  estimates  of  the  model  parameters.  These  two  input 
types  were  a  high-frequency  pitch  doublet  and  band-limited  white  noise. 

The  research  indicated  the  conclusion:  While  the  method  is  not  as  complex  as 
other  parameter  estimation  methods,  this  research  shows  that  linear  regression  can  be 
used  successfully  in  determining  longitudinal  dimensional  stability  derivatives  of  a 
tactical  missile  in  flight  when  using  a  control  input  form  with  higher  frequency 
modulations,  such  as  band-limited  or  filtered  white  noise. 
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PARAMETER  ESTIMATION  OF  A  MISSILE  USING  LINEAR  REGRESSION 


1.  Introduction 

1.1  General 

Aerodynamic  force  and  moment  databases  for  tactical  missile  trajectory 
simulations  are  usually  populated  with  an  evolving  mix  of  analytical  and  wind  tunnel 
derived  aerodynamic  predictions.  The  fidelity  of  these  simulations  can  be  improved  by 
incorporating  flight  test  derived  estimates  of  aerodynamic  characteristics.  Parameter 
estimation  allows  for  a  better  understanding  of  theoretical  predictions,  improves  wind 
tunnel  databases,  aides  in  the  development  of  flight  control  systems  and  provides  more 
accurate  representations  of  the  missile  in  all  flight  regimes. 

The  work  presented  here  is  concerned  with  the  parameter  estimation  of  a  tactical  missile 
utilizing  the  linear  regression,  least  squares  algorithm.  While  linear  regression  is  widely 
understood  in  the  scientific  and  mathematical  world,  little  work  has  been  done  in 
applying  the  least  squares  method  to  predict  aerodynamic  parameters  of  tactical  missiles. 

1.2  Background 

Several  methods  are  available  for  aerodynamic  system  identification  from  flight 
test  data.  The  primary  problem  is  determining  which  methods  are  most  applicable  to 
tactical  missiles.  Missile  flight  histories  are  characterized  by  rapid  accelerations,  rapidly 
changing  weight,  inertia  and  thrust  characteristics,  and  often,  short  flight  times.  These 
characteristics  make  accurate  parameter  estimation  of  missile  aerodynamics  more 
challenging  than  for  aircraft. 
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The  concept  of  aerodynamic  parameter  determination  from  flight  test  data  has 
been  successfully  applied  to  aircraft.  An  aircraft  used  in  flight  test  applications  can  be 
well  equipped  with  complex  instrumentation  and  sensors,  which  can  be  utilized 
repeatedly  for  numerous  tests.  Control  surface  inputs  can  be  designed  to  intentionally 
separate  complex  aerodynamic  effects  or  excite  aircraft  motion  while  at  a  particular  flight 
condition.  This  allows  for  the  usage  of  simplified  aerodynamic  models  linearized  about  a 
certain  flight  condition.  Flight  tests  can  also  be  repeated  if  necessary  when  sensor  or 
telemetry  problems  arise. 1 

In  contrast,  aerodynamic  parameter  determination  from  flight  test  data  for  a 
missile  flight  test  program  presents  a  more  difficult  problem.  Because  of  the  cost  and  the 
destructive  nature  of  a  missile,  relatively  few  flight  tests  can  be  afforded.  Flight  tests  that 
are  scheduled  usually  have  other  objectives  in  addition  to  detennining  the  aerodynamic 
parameters  of  the  missile.  Therefore,  the  opportunity  to  set  up  the  missile  control 
surfaces  to  determine  the  aerodynamic  parameters  might  not  be  available.  Most  missile 
flights  involve  large  and  rapid  variations  in  flight  condition.  This  minimizes  the  flight 
regime  in  which  a  simplified  aerodynamic  model  will  be  accurate  for  the  missile  flight 
test. 1  Therefore,  multiple  aerodynamic  models  are  needed  to  encompass  the  missile 
flight  test  envelope.  This  would  significantly  increase  computational  time  required  for 
the  missile  parameter  estimation  technique. 

Several  methods  used  for  system  identification  include,  Linear  Regression  (or  Least 
Squares),  Maximum  Likelihood  and  Extended  Kalman  Filters.  The  Maximum 
Likelihood  estimator  method  has  been  successfully  used  to  detennine  aerodynamic 
parameters  for  over  40  years.  Most  of  this  work  has  been  done  using  aircraft  and  not 
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much  research  has  been  found  applying  this  method  to  tactical  missiles.  The  Maximum 
Likelihood  method  is  popular  because  it  can  generate  statistically  optimum  estimates  of 
constant  parameters.  This,  however,  puts  the  Maximum  Likelihood  method  at  a 
disadvantage  when  estimating  missile  aerodynamic  parameters  due  to  the  fact  that  during 
the  missile  flight  test,  large  deviations  from  nominal  flight  conditions  are  to  be  expected. 
Popular  programs  that  implement  the  Maximum  Likelihood  estimator  are  PEST,  MMLE, 
and  PARAIDE. 

The  Extended  Kalman  Filter,  EKF,  method  generates  statistically  minimum  variance 

2 

optimum  estimates  that  can  be  time -varying  and  is  also  considered  numerically  efficient. 
This  allows  the  EKF  method  to  be  applied  to  flight  test  environments  in  which  the  system 
parameters  are  not  held  constant.  In  this  method,  the  system  equations  are  linearized 
about  the  state-parameter  vector  in  order  to  determine  the  optimal  system  estimate.  The 
maximum  likelihood  method  avoids  the  linearization  constraints  of  the  EKF  method; 
however,  by  doing  so,  it  introduces  excessive  computational  burdens  and  numerical 
difficulties. 

Linear  regression  has  been  used  recently  by  the  Air  Force  Research  Laboratory  for 
the  system  identification  of  on-line  reconfigurable  flight  control.  The  system 
identification  algorithm  studied  identified  rapid  changes  in  parameters  caused  by  failure 
and/or  damage.  Early  work  developed  relationships  from  flight  mechanics  to  regularize 
the  linear  regression  estimates.  Classical  linear  regression  was  augmented  to  include 
stochastic  constraints,  which  were  used  to  model  uncertainty  in  the  information  conveyed 
by  the  constraints.  In  the  interest  of  generality,  the  a  priori  estimates  of  the  stability  and 
control  derivatives  were  used  instead  as  regularization  constraints,  where  the  variance 
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represents  the  uncertainty  in  the  a  priori  estimates.  The  methodology  was  expanded  by 
additional  research  that  entailed  a  regularization  process  where  the  regressor  matrix  was 
subjected  to  a  singular  value  decomposition  that  led  to  a  reparametrization  of  the  original 
estimation  problem,  where  information  about  the  parameter  vector,  conveyed  by  the  data, 
was  extracted.  Then,  a  priori  information  was  used  to  regularize  the  estimation  problem 
and  obtain  an  estimate  of  the  original  parameter  vector.4  The  work  was  then  extended  to 
include  over-actuated  aircraft,  i.e.,  aircraft  with  distributed  control  effectors.  This 
development  was  demonstrated  on  the  F-16  VISTA  aircraft  for  both  the  lateral  and 
longitudinal  axes.  The  results  of  the  research  showed  that  even  during  periods  of  low 
excitation,  the  parameter  estimates  had  been  improved  by  utilizing  a  priori  infonnation. 

1.3  Research  Objective 

The  objective  of  this  research  is  to  detennine  the  applicability  of  using  the  Linear 
Regression  method  to  determine  the  aerodynamic  parameters  of  a  highly  maneuverable 
tactical  missile.  Specifically,  the  research  will  determine  the  longitudinal  aerodynamic 
parameters  and  compare  them  to  the  truth  data  to  detennine  the  viability  of  the  method 
for  use  as  a  parameter  identification  tool.  The  research  will  also  determine  the  best 
command  input  form  to  drive  the  6-Degree  of  Freedom  missile  simulation  that  will 
generate  good  estimates  from  the  Linear  Regression  method. 
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1.4  Approach  and  Scope 

The  work  here  presents  the  Linear  Regression  method  as  a  possible  tool  for 
parameter  identification  of  a  missile.  The  Linear  Regression  method  was  to  be  applied  to 
real  flight  test  data  from  a  missile  after  launch.  However,  data  from  real  flight  test 
experiments  was  not  available  for  this  research,  so  a  6  Degree-of-Freedom  missile 
simulation  in  MATLAB/SIMULINK  was  used  to  generate  a  data  set  representative  of 
data  collected  from  a  real  flight  test.  The  aerodynamic  data  for  the  6-DoF  missile  model 
was  based  on  a  supersonic,  tail  controlled  missile  with  similar  geometry  and  mass 
properties  of  an  AIM-9X  missile.  The  missile  aerodynamics  of  the  model  were  limited  to 
include  only  first-order  effects  due  to  limitations  in  the  semi-empirical  prediction  code, 
Missile  DATCOM,  utilized  for  this  research.  The  missile  model  does  not  include 
actuator  dynamics,  but  defines  different  combinations  of  individual  control  surface 
deflections  into  net  deflections  for  SP  (Roll  deflection  command),  5Q  (Pitch  deflection 
command)  and  SR  (Yaw  deflection  command). 

From  the  data  collected,  the  Linear  Regression  method  was  used  to  detennine  the 
three  longitudinal  dimensional  stability  derivatives,  Ma ,  M q  and  .  The  Linear 

Regression  method  was  validated  first  using  a  longitudinal  linear  model  of  the  system 
with  known  parameters.  These  parameters  were  then  estimated  using  the  Linear 
Regression  method  and  compared  to  the  known  parameters  for  validation. 

The  Linear  Regression  method  was  then  applied  to  the  parameters,  CL  (angle  of 
attack  in  degrees),  Cf  (pitch  rate  in  degrees/sec),  &  (Elevator  angle  in  degrees)  and 
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C{  (pitch  acceleration  in  degrees/sec")  obtained  from  the  6-DoF  simulation  data  to 
determine  estimates  of  the  longitudinal  stability  derivatives.  The  estimates  were 
compared  to  the  truth  data  from  the  aerodynamic  model  of  the  missile  simulation. 

The  Linear  Regression  was  applied  once  more  using  a  linear  model  of  the  system  that 
used  the  initial  estimates  obtained  from  the  simulation  data  as  the  new  truth  data.  Once 
again,  for  validation  purposes,  the  Linear  Regression  estimates  were  compared  to  the 
earlier  estimates  from  the  6-DoF  simulation. 

1.5  Relevance 

This  research  is  intended  to  investigate  the  feasibility  of  applying  Linear  Regression 
techniques  to  determine  aerodynamic  coefficients  of  a  missile  model.  It  is  not  intended 
to  compare  results  found  by  using  Linear  Regression  to  other  well-known  parameter 
identification  methods.  The  intent  of  this  research  is  to  provide  useful  findings  to  help 
improve  the  fidelity  of  aerodynamic  force  and  moment  databases  for  tactical  missile 
trajectory  simulations. 

1.6  Document  Overview 

The  document  starts  in  Chapter  2  with  an  introduction  to  Linear  Regression  and 
the  equations  and  their  derivations  that  were  used  in  this  research.  Chapter  3  provides 
the  background  needed  to  understand  the  missile  simulation  created  in  SIMULINK  that 
was  used  for  this  research.  It  includes  descriptions  of  the  missile  simulation,  missile 
aerodynamics  and  autopilot  that  make  up  the  6-degree-of-freedom  simulation.  Chapter  4 
combines  both  topics  from  Chapter  2  and  Chapter  3  to  produce  results  and  analysis  of  the 
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proposed  parameter  identification  method.  Chapter  5  goes  over  conclusions  found  from 
this  research  and  proposes  further  research  to  expand  on  what  was  learned. 
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2.  Linear  Regression  Technique 


2.1  Overview 

The  purpose  of  this  chapter  is  to  provide  a  review  of  the  theory  of  least  squares  or 
Linear  Regression,  as  is  described  from  Reference  5.  Least-squares  theory  has  become  a 
major  tool  for  parameter  estimation  of  experimental  data.  While  other  estimation 
methods  exist,  such  as  Maximum  Likelihood  and  Extended  Kalman  filters,  the  least- 
squares  method  continues  to  be  the  most  understood  among  engineers  and  scientists. 

This  is  due  to  the  fact  that  the  method  is  easier  to  understand  than  other  methods  and 
does  not  require  knowledge  of  mathematical  statistics.  The  method  also  may  provide 
solutions  in  cases  where  other  methods  have  failed.  Estimates  obtained  by  the  least- 
squares  method  have  optimal  statistical  properties;  they  are  consistent,  unbiased  and 
efficient. 

2.2  Least-Squares  Theory 

The  least-squares  technique  provides  a  mathematical  procedure  by  which  a  model 
can  achieve  a  best  fit  to  experimental  data  in  the  sense  of  minimum-error-squares. 
Suppose  there  is  a  variable  y  that  is  related  linearly  to  a  set  of  n  variables, 
x  =  (x1;x2,...,xn) ,  that  is 

y  ~  0\x\  +  $?x2  +  ...  +  0 nxn  (2.1) 

where  0  -  (6v62,...,6n)  is  a  set  of  constant  parameters.  9  are  unknown  and  we  wish  to 
estimate  their  values  by  observing  the  variables  y  and  x  at  different  times. 
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Assuming  that  a  sequence  of  m  observations  on  both  y  and  xhas  been  made  at  times 


,  we  can  now  relate  the  data  by  the  following  set  of  m  linear  equations: 

y(i)  =  0lxl(i)  +  O2x2(i)  + ...+ 0nxn(i),i  =  \,2,...,m  (2.2) 

This  equation  is  called  a  regression  function  and  9  are  the  regression  coefficients. 

This  system  of  equations  (2.2)  can  be  arranged  into  a  simple  matrix  form: 

y  =  xd  (2.3) 


where, 


>(  1) " 
y(  2) 
y(  3) 

x  = 

*i(!) 
A(2)  .. 

••  *„(1) 

••  x„(2) 

9  = 

A] 

92 

(2.4) 

_y(m)\ 

Xj  ( m )  . . 

..  xn(m)  \ 

A. 

There  are  several  different  instances  that  can  occur  when  determining  the  solution  of 
Equation  2.4  based  on  the  sizes  of  m  and  n. 

If  m  -  n ,  then  we  can  solve  for  6  uniquely  by 

0  =  X~ly  (2.5) 


provided  that  X  1 ,  the  inverse  of  the  square  matrix  X ,  exists.  6  denotes  the  estimate  of 
9. 

When  m  <  n  and  rank  ( A )  =  m  ,  there  are  an  infinite  number  of  solutions,  the  problem  is 
under-constrained.  In  general,  this  problem  has  an  infinite  number  of  solutions  9  which 
exactly  satisfy  y  -  X9  =  0  .  In  this  case  it  is  often  useful  to  find  the  unique  solution  for  9 
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which  minimizes  p|7 .  This  problem  is  referred  to  as  finding  a  minimum  nonn  solution  to 
an  underdetennined  system  of  linear  equations. 

In  the  most  usual  case,  m  >  n  and  rank(X)  =  n ,  there  is  no  solution  for  y  =  XO  .  The 
problem  is  also  referred  to  as  finding  a  least  squares  solution  to  an  over-detennined 

A  A 

system  of  linear  equations,  where  0  is  chosen  to  minimize  y  -  X  6  . 

The  best  approach  to  solving  the  over-determined  equation  (2.5)  would  be  to 

A 

determine  0  on  the  basis  of  least-error-squares  described  in  the  next  section. 

2.2.1  Least-Error-Squares 

For,  y  <£  R(X) 

Define  an  error  vector  e  =  (el,e2,...,em)1  and  let 

e  =  y  -  XO  (2.6) 

A 

6  will  be  chosen  in  such  a  way  to  minimize  the  criterion  J 

m 

J  =  fje?=eTe  (2.7) 

1=1 

To  carry  out  the  minimization,  J  is  expressed  as 

J  =  {y-XO)\y-XO) 

=  yTy  -  0TXTy  -  yT  XO  +  0T  XT  XO 

Differentiating  J  with  respect  to  0  and  equating  the  result  to  zero  to  detennine  the 

A 

conditions  that  the  estimate  0  minimizes  J,  provides  the  necessary  condition  for  the 
optimal  estimate 
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This  yields 


8J_ 

SO  eJe 


=  -2XTy  +  2XTX6  =  0 


XTXO  =  XTy  (2.8) 

A 

Solving  for  6 

e  =  {XTX)~lXTy  (2.9) 

Where  (XT Xyl  XT  is  commonly  referred  to  as  the  pseudo-inverse. 

This  result  is  known  as  the  least-squares  estimator  (LSE)  of  6  ,  also  called  ordinary  least 
squares. 

2.3  Statistical  Properties  of  Least-squares  Estimators 

In  this  section  the  statistical  qualities  of  the  least-squares  estimators  are 
examined.  Looking  at  Equation  (2.6),  where  the  vector  e  can  be  thought  of  as  the 
measurement  noise  and  or  modeling  error,  the  noise-disturbed  system  equation  is 

y  =  X6  +  e  (2.10) 

It  is  assumed  that  e  is  a  stationary  random  vector  with  zero  mean  value,  E \e\  =  0  and  that 
e  is  uncorrelated  with  y  andX ,  i.e.  white  noise.  Based  on  these  assumptions  about  e ,  it 
is  possible  to  determine  the  accuracy  of  the  parameter  estimates  given  by  equation  (2.9). 

A 

In  general,  0  is  a  random  variable.  The  accuracy  can  be  measured  by  a  number  of 
statistical  properties  such  as  bias,  error  covariance,  efficiency  and  consistency.  These 
tenns  will  be  defined  later  in  the  chapter. 
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A 

First  showing  that  9  is  unbiased,  meaning  that  E 


9 


=  9  .  Substituting  equation  (2.10) 


into  (2.9),  becomes 


9  =  9  +  (XTXylXTe 


Taking  the  expectation  on  both  sides  of  equation  (2.11)  and  applying  the 

A 

property  E [e]  =  0  ,  the  desired  result  is  obtained,  proving  that  9  is  unbiased. 


9 


=  E[9]  +  e][XtX) aXt  ]E[e]  =  9 


(2.11) 


(2.12) 


The  covariance  matrix  corresponding  to  the  estimate  error  9- 9  is 

x¥=E{(9-9)(9-9)t } 

=  E{[(Xt  xy1  XT e][(XT  xyl  XT e]T} 

=  ( XTXylXTE{eer}X(XTXy 1 
Define  the  covariance  matrix  of  the  error  vector  e  to  be 

R  =  E[eeT ]  (2.13) 

VF  is  reduce  to 

'F  =  (xT  xy1  xT  rx(xt  xy1  (2. 14) 

When  the  noise  e(i),  i  =  1,2,...,  are  identically  distributed  and  independent  with  zero 
mean  and  variance  cr2 ,  the  covariance  R  becomes 

R  =  E[eeT]  =  cr2I  (2.15) 
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This  implies  that  the  corresponding  LSE  0  is  a  minimum  variance  estimator.  6  is  called 
an  efficient  estimator. 


Lastly,  LSE  6  will  be  shown  as  a  consistent  estimator.  Rewriting  the  error  covariance 
matrix  'L  in  the  form  of 


w  =  a2(xTxy1=  — 

m 


1 


v1 


\m  j 


(2.16) 


where  R  =  a2 1  is  assumed  and  in  which  m  is  the  number  of  equations  in  the  vector 
equation  (2.10).  Assume  that  m)XT X\{  =  Y  ,  where  f  is  a  nonsingular 


constant  matrix.  Then 


=  lim,„„—  \-XTX 
m  ^  m  J 


=  0 


(2.17) 


Zero  error  covariance  means  that  6  =  0  at  m  — >  oo .  This  convergence  property  indicates 

A 

that  6  is  a  consistent  estimator. 

It  has  been  shown  in  this  section  that  the  LSE  in  the  presence  of  white  noise  is 
unbiased,  efficient,  and  consistent,  therefore  the  least  squares  technique  does  have  many 
advantages.  The  method  described  in  this  chapter  will  be  applied  to  results  obtained 
from  the  6-degree-of-Ereedom  missile  simulation  described  in  the  next  chapter,  Chapter 
3. 
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3.  Missile  Model  and  Simulation  Environment 


3.1  Overview 

This  chapter  presents  the  missile  simulation  that  will  be  used  to  generate 
simulation  data  that  the  Linear  Regression  technique,  described  in  Chapter  2,  will  employ 
for  parameter  identification.  The  missile  was  modeled  as  a  six-degree-of-freedom 
SIMULINK  model  within  a  1  vs.  0  simulation  developed  by  the  Air-to-Air  Weapons 
Branch,  Aircraft  and  Electronics  Division,  National  Air  and  Space  Intelligence  Center 
(ADNW/NASIC).  This  simulation  is  currently  used  to  detennine  missile  capabilities  and 
to  demonstrate  those  capabilities  in  enemy  engagement  scenarios.  This  section  discusses 
the  missile  data  and  the  simulation  obtained  from  ADNW/NASIC  and  changes  that  were 
made  to  facilitate  this  research. 

This  chapter  lays  out  the  background  needed  to  understand  the  Missile  simulation  created 
in  SIMULINK  that  was  used  for  this  research.  Section  3.1  describes  the  framework, 
components  and  environments  in  the  missile  engagement  simulation.  Section  3.2 
describes  the  aerodynamic  model  of  the  missile  first  by  defining  the  aerodynamic  angle 
and  control  surface  deflection  conventions.  The  aerodynamic  coefficients  and  forces  and 
moment  equations  are  then  presented  leading  up  to  finally  reach  the  Translational  and 
Rotational  Equations  of  motion  for  the  missile  airframe.  Section  3.2  ends  with  a 
description  of  the  type  of  missile  used  for  the  research  and  how  the  aerodynamic  data  was 
generated  for  the  model.  Section  3.3  describes  the  acceleration-  controlled  autopilot  used 
to  control  the  missile  model. 
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3.2  Missile  Simulation 


The  simulation  used  for  this  research  was  modified  from  an  existing  simulation 
developed  by  NASIC/ADNW  in  MATLAB/SIMULINK  that  was  based  on  the  missile 
trajectory  program,  TRAP  (TRajectory  Analysis  Program).  The  simulation  represented  a 
1  vs.  0  engagement  with  an  air-to-air  missile  fire. 

3.2.1  Simulation  Components 

The  simulation  is  setup  to  simulate  three  vehicles:  a  launch  aircraft,  a  missile  and 
a  target  aircraft.  It  is  built  around  a  detailed  fly-out  model  for  the  missile  with  simplified 
launch  aircraft  and  target  models.  The  launch  aircraft  is  modeled  as  a  pseudo  5-degree  of 
freedom,  or  modified  point-mass  model.  The  target  aircraft  is  a  3  degree  of  freedom  or 
point-mass  model  and  the  missile  is  a  6  degree  of  freedom  model.  All  three  models  in  the 
simulation  were  modeled  having  Flat-Earth  kinematics.  A  typical  engagement  is 
illustrated  below. 
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TYPICAL  ENGAGEMENT 


TARGET 

AIRCRAFT 


Figure  3-1:  Typical  Simulation  Engagement 


5  DoF  Launch  Aircraft 

The  launch  aircraft  is  modeled  as  a  modified  point-mass  with  angle-of-attack  and 
simplified  pitch  and  roll  dynamics  (also  referred  to  as  ‘pseudo  5-DOF’).  For  this 
research,  the  launch  aircraft  was  constrained  to  be  non-maneuvering  with  constant 
velocity  (maintained  straight  and  level  flight  path).  No  modifications  were  made  to  this 
model  for  the  research  described  in  this  paper. 

3  DoF  Target  Aircraft 

The  target  aircraft  is  modeled  as  a  simplified  point-mass.  The  target  maintained 
a  straight  and  level  flight  path  while  holding  constant  velocity.  This  setup  for  the  target 
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aircraft  was  sufficient  for  the  scope  of  this  research.  No  modifications  were  made  to  this 
model. 

6  DoF  Missile 

The  6  degree-of-freedom  missile  model  provided  by  ADNW  previously  had  not 
been  implemented  into  a  simulation  that  was  suitable  for  this  research.  The  6-DOF 
missile  model  replace  a  point-mass  missile  model  in  the  existing  1  vs.  0  simulation 
provided  by  ADNW.  The  original  DIME  Two  Axis  Gimbal  model  and  the  DIME  Simple 
Seeker  model  that  were  implemented  into  the  original  simulation  were  too  complex  for 
the  scope  of  this  research.  (DIME  stands  for  the  Air  Force  Research  Laboratories 
Munitions  Directorate  Dense  Inert  Metal  Explosives  Laboratory).  These  models  were 
replaced  with  a  much  less  complex  Momentum  Gimbal  model  and  a  Perfect  Seeker 
model,  respectfully. 

TELEMON  TRAP  Perfect  Seeker  and  TRAP  Momentum  Gimbal 

The  DIME  Simple  seeker  and  the  DIME  two-axis  gimbal  included  in  the  original 
simulation  were  overly  complicated  for  the  purpose  of  this  research.  These  models  were 
replaced  with  the  simple  TRAP  Perfect  Seeker  and  TRAP  Momentum  Gimbal.  For  the 
perfect  seeker,  the  seeker  axis  is  always  constrained  to  point  at  the  target  (subject  to 
gimbal  biases  and  limits)  and  the  seeker  line-of-sight  rate  commands  are  set  to  the  true 
line-of-sight  rate  resolved  into  vertical  and  horizontal  commands  in  the  body  axis  (the 
body  xz -plane  and  xy-plane).  The  seeker  line-of-sight  rate  commands  are  limited  to  the 
input  maximum  line-of-sight  tracking  rate.  The  momentum  gimbal  is  modeled  for  the 
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momentum-stabilized  seeker  platform.  The  seeker  x-axis  is  initially  constrained  to  point 
at  the  target  and  gimbal  angles  are  calculated  according  to  the  orientation  of  the  gimbal 
cages.  Gimbal  biases  are  input  in  pitch  and  yaw,  and  then  are  added  to  the  gimbal  angles. 
The  gimbal  angles  in  either  pitch  or  yaw  are  then  limited  to  the  input  maximum  gimbal 
angle.  The  resulting  gimbal  angles  are  used  to  calculate  a  new  seeker-reference-axes  to 
seeker-axes  transfonnation  matrix  and  hence  a  new  reference-axes  to  seeker-axes 
transfonnation  matrix.  This  matrix  is  used  to  determine  the  geometric  line-of-sight 
angular  tracking  errors  in  the  seeker.  For  the  perfect  seeker,  these  tracking  errors  are 
always  zero  when  the  gimbal  angles  are  not  limited  and  providing  the  gimbal  angle 
biases  are  zero. 

3.3  Missile  Aerodynamic  Model 

Now  that  the  overall  architecture  of  the  1  vs.  0  Missile  engagement  simulation  has 
been  introduced,  more  detail  is  going  to  be  presented  on  the  6-DoF  missile  model  used  in 
the  simulation. 

The  aerodynamic  methodology  applied  to  the  generic  missile  model  used  in  this 
research  is  based  upon  the  aerodynamic  methodology  from  the  NASIC/ADNA  software 
called,  the  Trajectory  Analysis  Program,  or  TRAP  for  short.  TRAP  has  been  widely  used 
to  evaluate  the  performance  of  aerodynamic  vehicles  in  l-vs-0  engagement  scenarios. 
Sections  3.3.1,  3.3.2  and  3.3.3  come  from  Reference  6. 
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3.3.1  Aerodynamic  Angle  Definitions 

Aerodynamic  characteristics  are  treated  as  a  function  of  total  angle-of-attack  (a 
and  aerodynamic  roll  angle  (</> ) .  The  relationship  between  these  angles  and  the  missile 


velocity  components  along  the  missile  body  x-axis  (MVELBX  =  u),  y-axis  (MVELBY  = 
v)  and  z-axis  (MVELBZ  =  w)  are  shown  in  Figure  3-2.  The  equations  relating  the  angles 
to  the  missile  velocity  components  are: 


</>' =  tan1-  (3.2) 

w 

Figure  3-2  also  shows  two  additional  angles,  these  are  PHIPRM  and  OFFPFII.  PFfTPRM 
is  the  aerodynamic  roll  angle  for  the  aerodynamic  characteristics  and  is  always  measured 
from  the  plane  containing  control  surfaces  1  and  3,  with  control  surface  1  to  leeward  for 
PHIPRM  =  0.0.  The  angle  OFFPHI  is  the  offset  angle  between  the  zero  body  roll  angle 
of  the  body  (in  this  case,  the  top  of  the  body  is  mid-way  between  controls  1  and  4)  and 
zero  aerodynamic  roll  angle  for  the  aerodynamic  characteristics.  OFFPHI  is  positive  if 
the  missile  must  be  rolled  in  the  positive  direction  (clockwise  looking  along  the  missile 
body  x-axis)  when  moving  from  zero  aerodynamic  roll  angle  (PHIPRM  =  0.0)  to  zero 
body  roll  angle  (PHIMSL  =  0.0).  Figure  3-2  shows  the  relationship  between  PHIPRM 
and  OFFPHI  for  both  a  ‘+’  and  an  ‘X”  configuration  missile  (OFFPHI  =  0.0  and  45.0 
degrees,  respectively). 
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Figure  3-2:  Aerodynamic  Angles 


3.3.2  Control  Surface  Deflection  Conventions 

The  aerodynamic  methodology  uses  the  concept  of  effective  pitch  and  yaw 
control  deflections  along  with  the  net  roll  control  deflection  and  the  net  squeeze  control 
deflection.  These  each  involve  different  combinations  of  individual  control  surface 
deflections.  The  conventions  used  by  the  aerodynamic  methodology  for  defining  control 
surface  deflections  are  described  in  the  sections  below  for  individual  control  surface 
deflections,  effective  pitch  control  deflection,  effective  yaw  control  deflection  and  roll 
control  deflection. 
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Individual  Control  Surface  Deflections 


The  convention  used  for  the  individual  control  surface  deflections  denotes  a 
positive  control  deflection  as  a  clockwise  rotation  of  the  fin,  from  looking  outboard  from 
the  missile,  see  Figure  3-3. 


PHIPRM  =  0.0  PHIPRM  =  45.0 

V  CONFIGURATION  V  CONFIGURATION 

Figure  3-3:  Individual  Control  Deflection  Convention 


Effective  Pitch  Control  Deflection 

The  effective  pitch  control  deflection  is  the  sum  of  the  individual  control  surface 
deflections  resolved  into  the  total  angle-of-attach  plane  according  to  the  equation: 

5Q  =  0.5(82  -  84)  cos  {PHIPRM)  +  0.5(fil  -  83)sm.(PHIPRM )  (3.3) 

where  DELQA  =  8Q  is  the  effective  pitch  control  deflection. 
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1 


PHIPRM  =  0.0  deg  PHIPRM  =  26.57  deg 

DELQA  =  20.0  deg  DELQA  =  22.36  deg 


4  1 


PHIPRM  =  45.0  deg 
DELQA  =  28.28  deg 


Figure  3-4:  Effective  Pitch  Control  Deflection 

Figure  3-4  shows  the  maximum  values  of  DELQA  when  the  missile  is  at  three  different 
aerodynamic  roll  angles  (PFIIPRM  =  0.0,  26.57,  45.0  degrees).  For  the  example  shown 
in  the  figure,  the  maximum  deflection  of  an  individual  control  is  20.0  degrees.  This  in 
turn  gives  the  maximum  values  of  DELQA  of  20.0,  22.36  and  28.28  degrees  at  the  three 
values  of  PHIPRM,  respectively.  The  maximum  values  follow  the  relationship: 

max  (A) 


ma  x(DELQA)  = 


cos  {PHIPRM) 


(3.4) 


where  max(A)  is  the  maximum  value  for  an  individual  control  surface  deflection.  For  the 
purposes  of  the  aerodynamic  tables,  it  is  inconvenient  to  have  different  maximum  values 
of  one  of  the  independent  variables  (DELQA)  associated  with  each  value  of  another  of 
the  independent  variables  (PHIPRM).  To  alleviate  this,  the  values  of  DELQA  are 
normalized  by  multiplying  the  actual  value  of  DELQA  by  the  cosine  of  PHIPRM.  The 
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illustration  in  Figure  3-4  shows  examples  of  individual  control  deflections  that  produce 
only  an  effective  pitch  control  deflection.  As  a  result,  the  magnitude  of  the  normalized 
value  will  never  exceed  the  maximum  magnitude  of  an  individual  control  surface 
deflection.  It  is  possible,  however,  to  obtain  normalized  values  of  the  effective  pitch 
control  deflection  that  are  greater  in  magnitude  than  the  maximum  magnitude  of  an 
individual  control  surface  deflection  and  exceed  the  maximum  value  in  the  data  table. 
This  occurs  when  the  demand  plane  is  not  aligned  with  the  total  angle-of-attack  plane.  In 
such  cases,  the  effective  pitch  control  data  will  be  extrapolated  to  the  required  value  of 
effective  pitch  control  deflection. 

Effective  Yaw  Control  Deflection 

The  effective  yaw  control  surface  deflection  is  the  sum  of  the  individual  control 
surface  deflections  resolved  nonnal  to  the  total  angle-of-attack  plane  according  to  the 
equation: 

SR  =  0.5(S3-S\)cos(PHlPRM)  +  0.5(S2-S4)sin(PHIPRM)  (3.5) 

Where  DELRA  =  SR  is  the  effective  yaw  control  deflection.  Figure  3-5  shows  the 
maximum  values  of  DELRA  when  the  missile  is  at  three  aerodynamic  roll  angles 
(PHIPRM  =  0.0,  26.57  and  45.0  degrees).  In  the  figure,  the  maximum  deflection  is 
assumed  to  be  20.0  degrees  which  gives  rise  to  the  maximum  values  of  DELRA  of  20.0, 
22.36  and  28.28  at  the  three  values  of  PHIPRM,  respectively. 

The  effective  yaw  control  deflection  is  never  used  as  an  independent  variable  in  the  data 
tables.  The  treatment  of  the  effective  yaw  control  deflection  is  handled  through  the  use 
of  aerodynamic  derivatives,  which  are  assumed  to  be  independent  of  the  magnitude  of  the 
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effective  yaw  control  deflection.  As  a  result,  there  is  no  need  to  normalize  the  values  of 
the  effective  yaw  control  surface  deflection.  This  use  of  aerodynamic  derivatives  can  be 
rationalized  because  the  magnitude  of  the  effective  yaw  control  deflection  will  generally 
be  small.  Large  effective  yaw  control  deflections  usually  occur  in  transient  situations  and 
are  of  short  duration. 


1  i 


PHIPRM  =  0.0  deg  PHIPRM  =  26.57  deg 

DELRA  =  20.0  deg  DELRA  =  22.36  deg 


1 


A 

PHIPRM  =  45.0  deg 
DELRA  =  28.28  deg 


Figure  3-5:  Effective  Yaw  Control  Deflection 
Roll  Control  Deflection 

The  net  roll  control  deflection  is  defined  by  the  equation: 

SP  =  0.25(S\  +  S2  +  S2>  +  84)  (3.6) 

Where  DELPA  =  8P  is  the  net  roll  control  deflection.  A  pure  roll  command  is  shown  in 
Figure  3-6  where  all  the  control  surfaces  are  deflected  by  equal  amounts  in  the  same 
direction  by  10.0  degrees. 
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A 

PHIPRM  =  0.0  deg 
DELPA  =  10.0  deg 


A 


PHIPRM  =  26.57  deg 
DELPA  =  10.0  deg 


A 

PHIPRM  =  45.0  deg 
DELPA  =  10.0  deg 


Figure  3-6:  Roll  Control  Deflection 

The  net  roll  control  deflection  is  never  used  as  an  in  dependent  variable  in  the  data  tables. 
The  treatment  of  the  net  roll  control  deflection  is  handled  through  the  use  of  aerodynamic 
derivatives,  which  are  assumed  to  be  independent  of  the  magnitude  of  the  net  roll  control 
deflection.  This  approach  is  rationalized  because  of  the  relatively  small  magnitudes  that 
are  usually  associated  with  the  net  roll  control  deflections. 

Squeeze  Control  Deflection 

There  is  a  combination  of  individual  control  surface  deflections  that  produce  no  net 
moment  and  no  lateral  force,  only  an  increase  in  the  axial  force.  This  is  called  the 
squeeze  control  deflection.  The  net  squeeze  control  deflection  is  defined  by  the  equation: 

SS  =  0.25{Sl-S2  +  S3-S4)  (3.7) 
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Where  DELSA  =  <55*  is  the  net  squeeze  control  deflection.  A  pure  squeeze  control  is 
illustrated  in  Figure  3-7,  where  the  maximum  deflection  for  each  individual  control 
surface  shown  is  10.0  degrees.  It  is  of  interest  for  control  system  design  to  keep  the  net 
squeeze  control  deflection  to  zero,  since  the  main  function  of  the  squeeze  control  is  to 
contribute  to  the  axial  force. 


PHIPRM  =  0.0  deg  PHIPRM  =  26.57  deg  PHIPRM  =  45.0  deg 

DELSA  =  10.0  deg  DELSA  =  10.0  deg  DELSA  =  10.0  deg 

Figure  3-7:  Squeeze  Control  Deflection 
Combinations  of  Control  Surface  Deflections 

The  relationship  between  the  individual  control  surface  deflections 
{dl,  52, 53, 5 A)  and  the  various  definitions  of  combined  control  surface  deflections 
(SP,SQ,SR,SS)  can  be  represented  by  collecting  the  equations  above  and  expressing 
them  in  the  matrix  form: 
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SP 

1 

1 

1 

1 

"£1 

SQ 

1 

2sin^ 

2  cos  (f> 

-  2  sin  ^ 

-2  cos  (j) 

S2 

SR 

~4 

-2  cos  </> 

2sin^ 

2  cos  (j> 

-2sin^ 

S3 

SS 

1 

-1 

1 

-1 

84 

Where  ^  is  defined  as  PHIPRM,  the  aerodynamic  roll  angle.  By  inverting  the  above 
matrix,  the  individual  control  surface  deflections  (SI,  52,  S3,  54)  can  be  determined  given 
any  set  of  combined  control  deflections  (SP,  SQ,  SR,  SS). 


S\ 

1  sin^ 

-COS  (j) 

1 

SP 

82 

1  COS  (j) 

sin^ 

-1 

SQ 

S3 

1  -  sin  (j) 

COS  (j) 

1 

SR 

S4_ 

1  -cos  (f) 

-siruj) 

-1 

SS 

(3.9) 


By  applying  the  constraint  of  zero  squeeze  control  deflection  to  the  above  equation,  the 
relationship  between  the  individual  control  deflections  and  three  axes  control  deflections 
is  determined: 


s\ 

1 

sin^ 

-COS  (j) 

52 

1 

COS (j) 

sin^ 

S3 

1 

-sin  (j) 

COS  (j) 

S4_ 

1 

-COS  (f) 

-sin^ 

SP 

SQ 

SR 


(3.10) 


3.3.3  Aerodynamic  Coefficients 

This  section  describes  the  aerodynamic  coefficients  that  are  included  in  the 
modeling  of  the  6-DOF  missile  used  in  this  research.  The  coefficients  described  in  this 
section  will  be  incorporated  into  the  forces  and  moments  equations  in  the  next  sections  of 
this  report. 
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Basic  Airframe  Aerodynamics 


There  are  six  aerodynamic  coefficients  that  describe  the  characteristics  of  the 
basic  airframe  without  any  control  deflection.  These  are: 

CA  Axial  force  coefficient 

CN a  Nonnal  force  coefficient,  measured  in  the  total  angle-of  -attack  plane 

CMREF '  ,  Pitch  moment  coefficient  about  the  reference  center-of-gravity  location, 

measured  in  the  total  angle-of-attack  plane. 

CT.  Side  force  coefficient,  measured  normal  to  the  total  angle-of-attack  plane. 

CLNREF^  Yaw  moment  coefficient  about  the  reference  center-of-gravity  location, 
measured  nonnal  to  the  total  angle-of-attack  plane. 

CLL  Roll  moment  coefficient 

Five  of  these  coefficients  are  a  function  of  the  aerodynamic  roll  angle,  total  angle-of- 
attack  and  Mach  number.  Note  that  the  subscript  a  denotes  that  the  coefficient  is 
referenced  to  the  total  angle-of-attack  plane. 

Airframe  Aerodynamics  with  Control  Effects 

There  are  six  aerodynamic  coefficients  and  12  static  aerodynamic  derivatives  that 
describe  the  incremental  aerodynamic  characteristics  due  to  control  deflection.  The  six 
aerodynamic  coefficients  are  associated  with  the  effective  pitch  control  deflection.  Six  of 
the  12  aerodynamic  derivatives  are  associated  with  the  roll  control  deflection  and  six  are 
associated  with  the  effective  yaw  control  deflection. 

Effective  Pitch  Control  Deflection 
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Because  the  missile  maneuvers  primarily  in  the  total  angle-of-attack  plane,  the 


effective  pitch  control  deflection  ( SQ  )  is  likely  to  be  any  value  up  to  the  maximum 


physical  limits  of  the  controls.  Therefore,  the  effect  of  the  effective  pitch  control 
deflections  is  treated  by  considering  the  actual  deflection  rather  than  using  an 
aerodynamic  derivative.  The  impact  of  effective  pitch  control  deflection  on  each  of  the 
six  basic  aerodynamic  coefficients  is  given  by: 

A  (CA\ 


a(«J, 


«  /  SQ 


A  [CMREFa  ) 


a  ’SQ 


a  (crj 


a  ’SQ 


Incremental  axial  force  coefficient  due  to  effective  pitch  control 
deflection,  SQ  . 

Incremental  normal  force  coefficient  due  to  effective  pitch  control 
deflection,  SQ  . 

Incremental  pitch  moment  coefficient  about  the  reference  center- 
of-gravity  location  due  to  effective  pitch  control  deflection,  SQ  . 

Incremental  side  force  coefficient  due  to  effective  pitch  control 
deflection  SQ  ,  measured  normal  to  the  total  angle-of-attack  plane. 


A(CLNREFa  )  Incremental  yaw  moment  coefficient  about  the  reference  center-of- 


gravity  location  due  to  effective  pitch  control  deflection  SQ  , 
measured  nonnal  to  the  total  angle-of-attack  plane. 


A  (CLL)^ 


Incremental  roll  moment  coefficient  due  to  effective  pitch  control 
deflection,  SQ  . 


The  use  of  the  prefix  A  indicates  that  this  is  an  incremental  value  in  the  coefficient  and 
the  suffix  SQ  indicates  that  it  is  due  to  a  specific  value  of  effective  pitch  control 
deflection,  SQ  .  Note  that,  by  definition,  the  total  angle-of-attack  can  only  be  positive 
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which  means,  for  a  statically  stable  missile  that  is  tail  controlled,  a  trim  condition  must 
always  be  associated  with  a  negative  value  of  SQ  .  These  coefficients  are  a  function  of 
the  effective  pitch  control  deflection,  total  angle-of-attack,  aerodynamic  roll  angle  and 
Mach  number. 


Roll  Control  Deflection 

Because  missiles  generally  have  small  roll  moments-of-inertia  compared  to  the 
aerodynamic  roll  moment  that  can  be  developed  by  the  control  surfaces,  usually  only 
small  control  surface  deflections  are  required  in  roll.  It  is  desired  then  to  use 
aerodynamic  derivatives  that  are  independent  of  the  magnitude  of  the  roll  control 
deflections  to  define  the  roll  control  effects.  The  effect  of  roll  control  deflection  (SP)  on 
the  aerodynamic  coefficients  is  given  by: 

Rate  of  change  of  the  axial  force  coefficient  with  the  square  of  the  roll 
control  deflection,  SP . 

Rate  of  change  of  the  normal  force  coefficient  with  the  square  of  the  roll 
control  deflection  SP ,  measured  in  the  total  angle-of-attack  plane. 

Rate  of  change  of  the  pitch  moment  coefficient  about  the  reference 

V  -  )sq 

center-of-gravity  with  the  square  of  the  roll  control  deflection  SP , 
measured  in  the  total  angle-of-attack  plane. 


f  dCA^ 

V  dP  j  gg 


f  dCN  ^ 


dCMREF . 

_ a_ 

r/P2 


dCY 

_ a_ 

dP 


J  SQ 


Rate  of  change  of  the  side  force  coefficient  with  the  roll  control 
deflection  SP ,  measured  nonnal  to  the  total  angle-of-attack. 
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f  dCLNREFa  i 


V 


dP 


JSQ 


Rate  of  change  of  the  yaw  moment  coefficient  about  the  reference 


center-of-gravity  with  the  roll  control  deflection  5P ,  measured  normal 
to  the  total  angle-of-attack  plane. 


Rate  of  change  of  the  roll  moment  coefficient  with  the  roll  control 
deflection  5P . 

The  use  of  the  suffix,  5Q  indicates  that  the  aerodynamic  derivative  applies  at  a  specific 
value  of  effective  pitch  control  deflection^  .  The  above  aerodynamic  derivatives  are  a 
function  of  the  effective  pitch  control  deflection,  total  angle-of-attack,  aerodynamic  roll 
angle  and  Mach  number. 


dCLL. 

dP 


JgQ 


Effective  Yaw  Control  Deflection 


When  the  maneuver  demand  plane  and  total  angle-of-attack  plane  are  not  aligned, 
the  result  of  resolving  the  individual  control  deflections  normal  to  the  total  angle-of- 
attack  plane  is  tenned  the  effective  yaw  control  deflection.  However,  because  the 
demand  plane  and  total  angle-of-attack  plane  are  generally  close  to  being  aligned,  the 
magnitude  of  the  effective  yaw  control  deflection  is  usually  very  small.  It  is  desired  then 
to  use  aerodynamic  derivatives  that  are  independent  of  the  magnitude  of  the  effective 
yaw  control  deflections  to  define  the  effective  yaw  control  effects.  The  effect  of  the 
effective  yaw  control  deflection  ( SR  )  on  the  aerodynamic  coefficients  is  given  by: 

UCA h 


dR- 


J  so 


Rate  of  change  of  the  axial  force  coefficient  with  the  square  of  the 
effective  yaw  control  deflection  SR . 
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dCN  , 


f  dCMREFa 


dR 2 


JSQ 


dCY 

_ a_ 

dR 


JSQ 


Rate  of  change  of  the  normal  force  coefficient  with  the  square  of 
the  effective  yaw  control  deflection  SR ,  measured  in  the  total 
angle-of-attack  plane. 

Rate  of  change  of  the  pitch  moment  coefficient  about  the  reference 
center-of-gravity  with  the  square  of  the  effective  control  deflection 
SR ,  measured  in  the  total  angle-of-attack  plane. 

Rate  of  change  of  the  side  force  coefficient  with  the  effective 
control  deflection  SR ,  measured  normal  to  the  total  angle-of-attack 
Plane. 


dCLNREF . 


dR 


JSQ 


Rate  of  change  of  the  yaw  moment  coefficient  about  the  reference 
center-of-gravity  with  the  effective  yaw  control  deflection  SR , 
measured  nonnal  to  the  total  angle-of-attack  plane. 


dCLNREF 

_ a_ 

dR 


JSQ 


Rate  of  change  of  the  roll  moment  coefficient  with  the  effective 
yaw  control  deflection,  SR . 


The  use  of  the  suffix  SQ  indicates  that  the  aerodynamic  derivative  applies  at  a  specific 
value  of  effective  pitch  control  deflection,  SQ  .  The  above  aerodynamic  derivatives  are  a 
function  of  the  effective  pitch  control  deflection,  total  angle-of-attack,  aerodynamic  roll 
angle  and  Mach  number. 
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The  coefficients  described  in  this  section  will  be  incorporated  into  the  forces  and 


moments  equations  in  the  next  section  of  this  report. 


3.3.4  Aerodynamic  Forces  and  Moments  Coefficient  Equations 

This  section  describes  the  steps  taken  in  the  missile  aerodynamic  model  to 
determine  the  aerodynamic  forces  and  moments  equations  used  to  derive  the  translational 
and  rotational  equations  of  motion  of  the  missile  airframe. 

First,  the  a  -  (j)  partial  aerodynamic  coefficients  are  combined  into  the  six 
primary  force  and  moment  coefficient  equations  in  the  wind  axis  system. 

The  total  aerodynamic  force  coefficient  terms,  including  both  longitudinal  and  lateral 
coefficients  are: 


CAtot  =  CA0  +  A(CAU  + 


fdCA\  2  JdCA^ 
(dP)  + 


V  dP  j  gQ 


V  dR  j  ^ 


m2 


(3.11) 


CY  =  CY  +  A(CY  .  )sn  + 

a  tot  a  0  a  ' 


f  dCY  ' 


V  dP  Jsq 


SP  + 


f dCY  ^ 


\  dR  Jsq 


SR 


(3.12) 


CN  ■  =  CN  .  +A(CA.)^  + 

a  tot  a  0  a  ' 


r  dCN  .  ^ 


m2  + 


f  dCN  ^ 


(SR)2  (3.13) 


The  total  aerodynamic  moment  coefficient  terms,  including  both  longitudinal  and  lateral 
coefficients  are: 


CLLtot  =  CLL0  +  A  (CLL)gQ  + 


dCLL 

dP 


SP  + 


Jsq 


dCLL 

dR 


SR 


(3.14) 


Jsq 


CMREF  =  CMREFq  +A(CMREF)S0  + 


>SQ 


dCMREF 

a 


\  r 

2 


dP2 


m+ 


Jsq 


dCMREF 

a 


dtf 


(Si)2  (3.15) 


Jsq 


44 


CLNREF  =  CLNREF  +A(CLNREF)jn  + 

n  a  o  V  n  '  cV 


>SQ 


r  dCLNREF ') 

a  1  SP+ 


dP 


JS2 


dCLNREF 

_ a_ 

dR 


Si  (3.16) 


JgQ 


Where  the  suffix  ‘0’  indicates  the  un-deflected  “basic  stability”  coefficient. 

Next,  the  six  primary  missile  aerodynamic  coefficients  described  in  Equations 
3.11  -  3.16  are  transfonned  from  the  a  -  (j)  Aeroballistic  coordinate  system  into  the 
Body  Axis  coordinate  system  with  the  use  of  the  coordinate  transformation  matrix: 


M  = 


cos  a  sin  a  sin  (j)  sin  a  cos  <f 
0  cos^  -sin^ 

-sin  a  cos  a  sin  (j)  cos  a  cos  <j> 


(3.17) 


Before  the  forces  and  moment  equations  can  be  developed,  the  dynamic  damping 
tenns  need  to  be  added  to  the  three  moment  coefficients  in  Equations  3.14-3.16.  These 
terms  consist  of  damping  derivative,  control  effectiveness  and  center  of  gravity,  c.g. 
adjustments.  Rewriting  Equations  3.14-3.16  with  these  tenns  gives: 


CLLtot  =  CLL0  +  A(CLL)rX)  +C,^-  + 


2V 


dCLL 

dP 


8P  + 


JSQ 


dCLL 

dR 


SR 


(3.18) 


JSQ 


CMREF  ,  =  CMREF  ,  +  A(CMREF  .)m+C  —  + ^  +  ... 

"  a  0  v  a  JdQ  mq  -- 


qj_  ,  CNaK,ref  ~  *  J 

2V 


l 


dCMREF  ,i  ,  r 
tt  1  (SP)2  + 


dP- 


JgQ 


dCMREF  ,i 

a  1  (SR)2 


v 


dR 


J  5Q 


(3.19) 
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rl  CY.lx  f—x  ) 

CLNREF  =  CLNREF.  +A (CLNREF.)  m+C  —  +  gV  cgM — c-^+. 

a  tot  a  0  v  a '  nr  J 


(3.20) 


r  dCLNREF  \ 

a  1  SP+ 


v 


dP 


r dCLN RE F\ 

a  1  SR 


JSQ 


v 


dR 


JSQ 


Where, 


Cj  Variation  of  rolling  moment  coefficient  with  rate  of  change  of  roll  rate 
Cm  Variation  of  pitching  moment  coefficient  with  pitch  rate 
Cn  Variation  of  yawing  moment  coefficient  with  rate  of  change  of  yaw  rate 
p  Roll  rate  of  the  missile 

q  Pitch  rate  of  the  missile 

r  Yaw  rate  of  the  missile 

/  Missile  aerodynamic  reference  length 

V  Airspeed 

xcg  ,-ef  Reference  Center-of-gravity  of  nose  where  the  Aerodynamic  tables  were 

generated  from. 

Center-of-gravity  of  missile 


3.3.5  Accounting  for  First  Order  Effects 

Before  continuing  on  to  the  forces  and  moments  equations,  the  forces  and 
moment  coefficient  terms  must  be  simplified  to  include  only  first  order  effects.  The 
aerodynamic  data  used  in  the  missile  model  does  not  include  second-order  effects 
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because.  Missile  DATCOM,  the  semi-empirical  prediction  software  used  to  create  the 
aerodynamic  data  is  not  capable  of  generating  any  second-order  control  effects. 

The  first-order  aerodynamic  characteristics  are  classed  as  those  in  which  the 
control  surface  deflection  has  a  primary  effect  on  the  particular  aerodynamic  coefficient. 
The  first-order  effects  are  the  effect  of  any  of  the  three  control  deflections  {5Q  ,8P ,8R) 
on  the  axial  force  coefficient,  the  effect  of  5Q  on  the  normal  force  and  pitch  moment 
coefficients,  the  effect  of  §P  on  the  roll  moment  coefficient  and  the  effect  of  5R  on  the 
side  force  and  yaw  moment  coefficients.6 
The  first-order  terms  are  as  follows: 

Basic  stability  coefficients: 

CA,CN  ,  CMREF  CY.,  CLNREF . ,  CLL 

a  a  a  a 

Incremental  coefficients: 


A(C4)*,  A(CA  ,  )*,  A  {CMREF a  )* 


Aerodynamic  derivatives: 


dCA 
~dP 2 


JsQ 


ClCLL 
v  dP 


JgQ 


dCA 

v  dR  j  gQ 


dCY . 

_ a_ 

dR 


JgQ 


dCLNREF , 


dR 


JSQ 


Equations  3.11-3.13  and  3.18-3.20  are  re-written  to  include  only  the  first-order  terms: 
For  the  total  force  coefficient  terms: 


CAtot=CA0+A(CA)SQ  + 


dCA 
dP 2 


(8P)2  + 


JgQ 


dCA 

dR2 


m2 


J  SQ 


(3.21) 


f dCY . ) 

CY  =  CY  .  +  - SR 


a  0 


V 


dR 


JgQ 


(3.22) 
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CN  .  =CN  .  +A (CN  .)m 

a  tot  a  0  a  ' 


(3.23) 


For  the  total  moment  coefficient  terms: 


CLL,,,=CLL„+ClrJ^  + 


f  dCLO 


V  dP  j  gQ 


8P 


(3.24) 


CMREF,  =  CMREF  .  +  A(CMREF  +  C  -^-  + 

n  tnt  (X  0  ^  ry  '  0(J  Wl„  r\xr 


ql  ,  CNa{xcg^ref-xcg) 


2V 


(3.25) 


CLNREF  =  CLNREF  +C  —  + 

a  tot  a  0  "i-  2Y 


rl  CY  (x . ,  -x„. 


eg, ref  eg) 


/ 


+ 


d  CLNREF 


dR 


8R  (3.26) 


JdQ 


3.3.6  Equations  of  Motion 

The  equations  of  motion  presented  here  are  for  a  missile  over  a  flat  Earth.  Now 
that  the  force  and  moment  coefficients  have  been  defined  with  first-order  terms,  the 
forces  and  moments  equations  now  can  be  presented. 

The  force  equations  are: 

Fx  =  - qSCA, ot  +  THR  UST  (3.27) 


(3.28) 

(3.29) 

The  moment  equations  are: 


F,  =  qSCY , 

1  1  a  tot 


F3  =  -qSCN 


Mx  =  qSdCLL 


(3.30) 
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M2  =  qSdCMREF 


(3.31) 


M3  =  qSdCLNREF . 


(3.32) 


Where, 

q  Dynamic  pressure 

S  Reference  Area 

d  Missile  Diameter 

Now  the  traslational  and  rotational  equation  of  motion  for  the  missile  model  can 
be  presented.  The  translational  degrees  of  freedom,  represented  by  the  velocity 
components  u  ,  vand  w,  are  solved  by  Newton’s  equation;  and  the  rotational  DoF, 
which  are  expressed  in  body  rates  p  ,  q  and  r ,  are  governed  by  Euler’s  equation.7 
Newton’s  second  law  with  respect  to  Earth,  as  the  inertial  reference  frame,  states  that  the 
time  rate  of  change  of  linear  momentum  equals  the  externally  applied  forces.  Newton’s 
Equations  for  the  translational  DoF  are: 


du  F. 

—  =  rv-qw-\ — -  +  t  g 
dt  m  13 


dv  F1 

—  =  pw  -ru  ^ — -  +  t  g 
dt  m  23 


du 


F 


—  =  qu-pv  +  —  +  tg 
dt  m 


(3.33) 

(3.34) 

(3.35) 


Where, 

t  Element  in  the  direction  cosine  matrix 


g  Acceleration  due  to  gravity 
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Euler’s  law  states  that  the  time  rate  of  change  of  angular  momentum  equals  the  externally 
applied  moments.7  Euler’s  Equations  for  the  rotational  DoF  are: 


J  =  V'fc-^zzV  +  ^l) 

(3.36) 

^  =  v'((4  -0/"'  +  we) 

(3.37) 

^  =  IzZ~l((lXX-Iyy)p<l  +  m  3) 

(3.38) 

3.3.7  Aerodynamic  Data  Generation 

The  aerodynamic  data  for  the  6-DoF  missile  model  was  based  on  a  supersonic, 
tail  controlled  missile  with  similar  geometry  and  mass  properties  of  an  AIM-9X  missile. 
The  data  was  generated  using  the  semi-empirical  prediction  code,  Missile  DATCOM,  or 
MISDATA. 

This  section  discusses  the  use  of  the  aerodynamic  prediction  method,  Missile 
DATCOM,  or  MISDAT,  to  generate  the  aerodynamic  table  data  required  for  the  6 
degree-of- freedom  SIMULINK  missile  model  as  presented  by  Reference  6.  The  fidelity 
associated  with  semi-empirical  prediction  methods  is  such  that  in  most  instances,  the 
aerodynamic  characteristics  will  repeat  at  every  90.0  degrees  in  aerodynamic  roll  angle. 
This  considerably  reduced  the  amount  of  data  that  was  required  to  represent  the  airframe 
aerodynamic  characteristics. 
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Basic  Airframe  Aerodynamics 


The  basic  airframe  aerodynamics  consists  of  data  on  the  missile  configuration  without 
the  controls  being  deflected.  Since  the  missile  is  modeled  in  the  cruciform  shape 
(PHIPRM  =  0.0)  and: 

•  Is  assumed  to  have  a  small  degree  of  asymmetry 

•  Does  not  operate  at  a  high  angle-of-attack, 

The  aerodynamic  characteristics  are  well  enough  behaved  to  be  repeated  every  90.0 
degrees  in  aerodynamic  roll  angle.  Since  this  was  the  case,  the  data  was  obtained 
through  a  range  of  45.0  degrees  in  aerodynamic  roll  angle,  which  captured  a  complete 
half  cycle  of  the  data.  The  number  of  intennediate  roll  angles  depends  on  the  behavior  of 
the  aerodynamic  characteristics.  The  use  of  a  small  number  of  intermediate  values  can  be 
accepted  where  the  data  are  well  behaved  when  combined  with  the  use  of  Hermitian 
interpolation.  The  smallest  increment  used  in  the  data  set  was  +  5.0  degrees.  A 
complete  listing  of  all  the  basic  airframe  aerodynamic  variables  and  the  associated  value 
ranges  they  were  collected  at  is  listed  in  Table  1 . 

Effect  of  Controls 

For  basic  missile  configurations  that  only  exhibit  a  small  degree  of  asymmetry, 
the  incremental  effect  of  the  control  deflections  is  very  small.  Therefore,  for  a  cruciform 
missile  with  little  asymmetries,  the  control  effectiveness  is  provided  over  a  45.0  degree 
range  of  aerodynamic  roll  angle. 

Since  it  is  unlikely  that  semi-empirical  prediction  methods  are  capable  of  generating  any 
of  the  second-order  effects  of  the  controls,  the  data  set  used  in  this  research  was  limited 
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to  first-order  effects  only.  The  first  order  effects  are  those  in  which  the  change  in  the 
aerodynamic  coefficient  is  in  the  same  plane  in  which  the  control  deflection  occurs.  The 
first-order  effects  due  to  effective  pitch  control  deflection  are  actual  increments  in  the 
values  of  the  aerodynamic  coefficients.  For  the  roll  control  deflection  and  effective  yaw 
control  deflection,  the  first  order  effects  are  represented  by  aerodynamic  derivatives. 

A  complete  listing  of  all  the  effects  of  control  variables  and  the  associated  value  ranges 
they  were  collected  at  is  listed  in  Table  1. 
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Table  3-1:  Aerodynamic  and  Control  Derivatives 


Variable 

Description 

Function  of: 

Mach 

Phi  (deg) 

Alpha 

(deg) 

DelQ 

Xcg 

ACAO 

Undeflected  Axial  Force 

ACNO 

Undeflected  Normal  Force 

ACYO 

Undeflected  Side  Force 

(j)  ,a  ,M 

0.4  to  5.0 

0  to  45 

0  to  25 

N/A 

N/A 

ACMRFO 

Undeflected  Pitching  Moment 

ACNRFO 

Undeflected  Y aw  Moment 

ACLLO 

Undeflected  Roll  Moment 

Variable 

Description 

Function  of: 

Mach 

Phi  (deg) 

Alpha 

(deg) 

DelQ 

Xcg 

ADLCAQ 

Incremental  Axial  Force 

ADLCAP 

Axial  Force  Deriv.  (wrt  da2) 

ADLCAR 

Axial  Force  Deriv.  (wrt  dr2) 

ADLCNQ 

Incremental  Normal  Force 

ADLCNP 

Normal  Force  Deriv.  (wrt  da2) 

ADLCNR 

Normal  Force  Deriv.  (wrt  d,2) 

ADLCYQ 

Incremental  Side  Force 

ADLCYP 

Side  Force  Deriv.  (wrt  da) 

5e,(f  ,a  ,M 

0.4  to 
5.0 

[0  26.5  45.0] 

0  to  25 

-20  to  20 

N/A 

ADLCYR 

Side  Force  Deriv.  (wrt  dr) 

ADCMRQ 

Incremental  Pitching  Moment 

ADCMRP 

Pitch  Moment  Deriv.  (wrt  da2) 

ADCMRR 

Pitch  Moment  Deriv.  (wrt  dr2) 

ADCNRQ 

Incremental  Y aw  Moment 

ADCNRP 

Y aw  Moment  Deriv.  (wrt  da) 

ADCNRR 

Y aw  Moment  Deriv.  (wrt  dr) 

ADCLLQ 

Incremental  Roll  Moment 

ADCLLP 

Roll  Moment  Deriv.  (wrt  da) 

ADCLLR 

Roll  Moment  Deriv.  (wrt  dr) 

Variable 

Description 

Function 

of: 

Mach 

Phi’ 

(deg) 

Alpha 

(deg) 

DelQ 

Xcg 

CMQ 

Pitch  Damping  (wrt  pitch 
rate) 

M,  XCG 

0.4  to  5.0 

N/A 

N/A 

N/A 

[1.52  1.58  1.65] 

CMAD 

Pitch  Damping  (wrt 
alphadot) 

M,  XCG 

0.4  to  5.0 

N/A 

N/A 

N/A 

[1.52  1.58  1.65] 

CLP 

Roll  Damping 

M 

0.4  to  5.0 

N/A 

N/A 

N/A 

N/A 
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3.4  Missile  Autopilot 

The  original  missile  model  did  not  have  an  acceleration  autopilot  implemented;  in 
order  to  control  the  missile  during  simulation,  pitch  plane  and  yaw  plane  acceleration 
autopilots  were  added.  The  acceleration  autopilots  were  taken  from  Zipfel,  Reference  7. 
Figure  3-8,  obtained  from  Reference  7,  shows  a  general  diagram  of  the  acceleration 
autopilot  used  in  the  missile  simulation.  Because  this  is  a  missile  simulation,  the  yaw 
plane  is  implemented  in  the  same  way. 

The  acceleration  tracking  loops  for  both  autopilots  were  designed  using  modern  pole 
placement  techniques.  For  best  performance,  proportional  and  integral  (PI)  techniques 
were  applied.  Proportional  control  was  used  for  quick  response  and  integral  control  was 
used  for  zeroing  the  steady-state  errors.  To  improve  performance,  an  inner  rate  loop  was 
added  for  stability  augmentation.  The  three  feedback  gains  in  the  autopilots  were  solved 
to  satisfy  the  specified  closed-loop  response. 


Figure  3-8:  General  Acceleration  PI  Autopilot 

The  linear  time-variant  plant  is: 

x  =  F  (t)x  +  g(t)u  (3.39) 
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From  figure  3-8,  the  proportional  and  integral  feed-forward  branches  with  their 
respective  gains,  GPand  GI  are  easily  visible.  The  major  feedback  loop  is  via  the  rate 

and  acceleration  gains  c  =  [k2  &,].  A  second  acceleration  feedback  with  unit  gain  lu  is 

wrapped  around  the  outside  to  improve  performance.  The  relationship  for  the  control 
variable  can  be  derived  from  the  figure. 

u  =  -cx  +  Gj  J  (v  -  h,x)dt  +  Gp ( v  -h,x)  (3 .40) 

The  states  are  then  augmented  by  introducing  the  scalar  auxiliary  variable;  in  its  state 
equation  form: 


£  =  v-Iux  (3.41) 

Substituting  u  into  the  open-loop  system  yields  the  closed-loop  system,  augmented  by  the 
auxiliary  variable^ : 

^1  —  rr( /■  -1-  Cl  li  \  Cl  rr~lr  v~l  \  Cl  cr~ 1 

(3.42) 


X 

F  -  g(c  +  GPF)  G,g 

X 

-L 

\GPg] 

v 

4 

i 

i 

Sr- 1 

| 

o 

1 _ 

i 

V 

l 

Which  can  be  abbreviated  as: 


x  =  F  (t)x  +  g  (t)v  (3.43) 

The  eigenvalues  of  this  closed-loop  fundamental  matrix  F  must  be  equal  to  the  desired 
closed-loop  poles.  The  condition  for  pole  placement  is: 


det 


4^  -  F  +  g(c  +  GP  F )  -  Gjg 
h =  s 


=  fl(‘s'”  A) 


(3.44) 


i=i 


Given  the  linearized  longitudinal  plant  equation: 
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q 

M  ^ 
q  N 

a 

q 

+ 

1 

1 _ 

a 

I 

i 

l _ 

a 

L  0  J 
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(3.45) 


where  a  is  the  normal  acceleration,  and: 


N„ 


M 


M 


Sq 


qs  r 

m  N° 

(3.46) 

qSd2 

2 12V  m“ 

(3.47) 

qSd 

1 

1  2 

(3.48) 

qSd 

T 

1  1 

(3.49) 
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For  the  condition  of  pole  placement,  the  corresponding  equations  are: 

x  =  \q  a\ 


F  = 


N„ 


M , 


Na 

Nr 


V 


(3.50) 


(3.51) 


g  = 


0 


(3.52) 


u  =  8 


(3.53) 


A,=[0  l] 


(3.54) 
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c  =  \k2  k] 


(3.55) 


Figure  3-9  Acceleration  Autopilot  for  Missile 


Figure  3-9  shows  in  detail,  the  variables  described  above  added  to  Figure  3-8  (from 
Reference  7).  This  schematic  was  directly  implemented  into  the  model  simulation 
SIMULINK  code.  The  gains  kx ,  k2  and  Gt  are  calculated  from  the  pole  placement 
condition  equation  (3-11).  Evaluating  the  left-hand  determinant  and  equating  terms  of 
equal  power  yields  the  three  gains: 


Gj  = 


co2  p 
NaMs 


k2  = 


M 


8  V 


2  (co  +  p  +  Mq-^- 


k ,  = 


1  NaMs 


MnN„  ,  MSN. 


a>~  +  2  £cop  +  Ma  +  ■  ,v2 


8ly  a 


V 


(3.56) 

(3.57) 

(3.58) 
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The  autopilot  was  implemented  into  the  simulation  to  calculate  these  gains  on-line  with 
changing  airframe  conditions.  The  system  dynamic  characteristics,  natural  frequency, 
damping  and  pole  location  values  chosen  were  evaluated  over  the  missile  flight  envelope 
to  ensure  stability  and  desired  system  behavior.  The  solution  for  the  position  feed¬ 
forward  gain  Gp  is  not  accessible  by  the  pole  placement  technique  and  was  determined 
based  on  root  locus  analysis  of  the  system. 

To  summarize,  this  chapter  laid  the  background  needed  to  understand  the  Missile 
simulation  created  in  SIMULINK  that  was  used  for  this  research.  Chapter  4  will  now 
utilize  outputs  generated  from  running  the  simulation  along  with  the  linear  regression 
methods  described  in  Chapter  2  to  analyze  linear  regression  as  a  tool  for  parameter  ID  of 
a  missile. 
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4.  Results  and  Analysis 


4.1  Aerodynamics  of  the  Experiment 

This  research  focused  on  the  longitudinal  equations  of  motion,  primarily  the  pitch 
equations  of  motion.  The  uncoupled  pitch  dynamics  consist  of  the  pitching  moment 
equation  and  the  normal  force  equation.  The  longitudinal  angular  equation  for  pitch 

acceleration,  or  <7  is: 

q  =  Maa  +  Mqq  +  M&Se  (4.1) 

and  the  corresponding  normal  force  equation  is: 

a  =  Zaa  +  Zqq  +  Z&5e  (4.2) 

Combining  equations  4. 1  and  4.2  into  a  linear  state-space  model  format  gives: 


The  remainder  of  the  research  concentrates  on  solving  for  Ma ,  M  and  M&  from 

Equation  4.3  using  the  method  of  Linear  Regression  and  comparing  those  results  to  “truth 
data”  obtained  from  the  simulation. 

4.2  Linear  Regression  Method  Validation 

Before  the  experimental  research  began,  the  Linear  regression  method  described 
in  Chapter  2  was  validated  using  the  linear  state-space  model,  Equation  4.3,  with  a 
simple  set  of  longitudinal  stability  and  control  derivatives.  This  data  was  taken  from 
Reference  8,  for  a  small,  single  jet  engine,  military  training  airplane. 
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For  an  aircraft,  the  equations  for  the  longitudinal  dimensional  stability  derivatives 


(from  Reference  8)  are: 


Z  = 


~qS(C,  +CD) 


m 


-qScC, 

Z  = _ : _ A 

q  2  mU 


7  _-qSCLSe 
m 


= 


qScCm 


„  qs ~c2c 

M  = - 5 

9  21  U 


(4.4) 

(4.5) 

(4.6) 

(4.7) 

(4.8) 


M&  =  qScC mSe  (4.9) 

yy 

For  the  example  aircraft,  the  following  values  were  used  to  solve  the  above  equations. 
Flight  condition: 

25,000  ft.  altitude,  Mach  0.6,  Cruise 
Reference  Geometry: 

S  (wing  area)  =  136  ft2,  c  (root  chord)  =  4.4  ft 
Flight  Condition  Data: 

U  (True  Airspeed,  TAS)  =  610  ft/sec 
q  (Dynamic  Pressure)  =198  lbs/ft" 

Mass  Data: 

W  (weight)  =  4,000  lb 
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Iyy  =  4,800  SlUg-ft2 


Steady  State  Coefficients: 

C,  =0.149 
CD  =  0.0220 

Longitudinal  Coefficients  and  Stability  Derivatives  (Stability  Axis,  Dimensionless): 
C,  =  5.5 

CLq  =10.0 
=  -0.24 
C,„9  =-17.7 

Longitudinal  Control  and  Hinge  Moment  Derivatives  (Stability  Axis,  1/rad) 

C,  =0.38 

Cm  =-0.88 

mSe 


Solving  for  Equations  4.4  through  4.9  and  putting  their  values  into  Equation  4.3  gives, 


-37.17  -0.2979 

■290.82  -2895.5 


2.558 

26.66 


Se 


(4.10) 


The  linear  state-space  model  was  then  simulated  in  MATLAB  using  the  command 
“LSIM”  using  a  sine  wave  with  amplitude  of  100  as  the  control  elevator  input  for  5 
seconds.  The  resulting  time  histories  for  Angle  of  Attack,  Pitch  rate,  Elevator  deflection 
and  Pitch  acceleration  for  the  simulation  are  shown  below  in  Figure  4-1. 


61 


Alpha 


Elevator  Deflection 


Pitch  Rate 


Pitch  Acceleration 


Figure  4-1  Time  Histories  of  Linear  model  parameters 

To  validate  the  linear  regression  method  described  in  Chapter  2,  the  time  history  values 
of  a  (angle  of  attack  in  degrees),  q  (pitch  rate  in  degrees/sec),  Se  (Elevator  angle  in 
degrees)  and  q  (pitch  acceleration  in  degrees/sec“)  were  obtained  from  the  linear 
simulation  to  estimate  the  derivatives,  Ma ,  M  q  and  using  Equation  4.1. 

q  =  Maa  +  Mqq  + 

Using  Equation  2.9  from  Chapter  2  for  the  least-squares  estimator 

0  =  (XTX)-lXTy 

and  defining  X  and  Y  as: 

X  =  \a  q  Se\  (4.1 1) 
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(4.12) 


Y  -  q 

A 

The  values  obtained  for  6  using  Equation  4.1 1,  4.12  and  2.9  are  compared  to  the  original 
values  calculated  from  Equations  4.4  through  4.9  in  Table  4-1  below. 


Table  4-1:  Comparison  of  Longitudinal  Dimensional  Stability  Derivatives 


0  Values 

“Truth  Data” 
(Equations  4.4  -  4.9) 

Ma 

-290.8 

-290.82 

Mq 

-2895.5 

-2895.5 

M5e 

-26.70 

-26.66 

From  Table  4-1,  it  is  shown  that  the  linear  regression  method  produces  results  that 
compare  well  with  “Truth  Data”.  Slight  variations  in  the  numbers  are  due  to  rounding  in 
MATLAB.  This  validates  the  method  for  the  scope  of  this  research. 

4.3  Exercising  the  Simulation 

The  6-DOF  missile  simulation  described  in  Chapter  3  was  exercised  to  create  a 
series  of  batch  runs  containing  several  different  missile  fly-out  scenarios  terminating 
when  the  missile  reached  the  target.  After  reviewing  the  accumulated  data  from  the 
simulation,  it  was  found  that  the  autopilot  of  the  missile  model  was  not  commanding  a 
suitable  input  for  excitation  of  missile  maneuvers  used  for  parameter  identification. 

4.3.1  Control  Signal  Input  Form 

Determining  the  best  input  form  to  successfully  excite  the  missile  airframe  for 
parameter  identification  study  is  one  of  the  most  critical  parts  of  the  design  of  the 
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experiment.  Many  different  fonns  of  control  inputs  have  been  used.  There  are  usually 
four  common  kinds: 

•  Transient  (Step,  Impulse  etc.) 

•  Sinusoids  of  various  frequencies  and  amplitudes  (Doublets) 

•  Stochastic  signals 

•  Pseudorandom  (white)  noise 

It  is  well  known  that  the  shape  of  an  input  signal  could  influence  the  accuracy  of 
estimated  parameters  from  dynamic  flight  measurements.9 

The  power  of  an  input  signal  should  be  distributed  uniformly  over  the  frequency  range 
covering  the  dominant  airframe  dynamics.  Transient  signals  such  as  a  step  input  or  an 
impulse  usually  fail  in  this  regard  due  to  their  short  duration  over  the  simulation.  All 
modes  of  the  system  might  not  be  excited  (lack  of  persistent  excitation)  and  the  regressor 
matrices  used  to  obtain  model  parameters  might  be  close  to  being  singular  leading  to 
numerical  problems.  In  Figure  4-1  taken  from  Reference  9,  power  spectral  densities  of 
three  inputs  are  presented.  The  step  input  excites  modes  of  lower  frequencies  only, 
which  makes  it  unsuitable  for  parameter  estimation.  The  power  spectral  density  of  a 
multi-step  input  is  a  relatively  wide  band  compared  to  the  doublet,  which  excites  only  a 
particular  band  at  a  higher  frequency.  By  changing  the  duration  of  the  doublet,  the  peak 
of  the  power  spectral  density  can  be  shifted  to  lower  or  higher  frequencies,  which  makes 
it  very  suitable  for  parameter  estimation. 

System  identification  using  pseudorandom  noise  signals  is  also  very  commonly 
used.  The  noise  is  usually  rich  in  all  frequencies  and  therefore  excites  all  modes  of  the 
plant  being  measured.  This  leads  to  good  estimates  of  the  model  parameters.10 
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INPUT  SIGNALS: 


Figure  4-2  Frequency  Domain  Comparison  of  Various  Input  Signals 

Klein,  Reference  9,  has  determined  other  basic  requirements  for  aircraft  parameter 
identification: 

•  The  input  form  should  be  selected  in  agreement  with  the  mathematical  model 
representing  the  aircraft  under  test.  For  example,  an  input  for  the  longitudinal  short- 
period  model  with  linear  stability  and  control  derivatives  should  not  cause  aircraft  motion 
where  the  assumptions  of  constant  airspeed  and  linear  aerodynamics  are  not  valid. 

•  The  input  form  should  also  be  within  the  bandwidth  of  the  actuator  driving  the 
control  surface  being  commanded  for  the  maneuver  of  the  missile.  The  missile  model 
used  for  this  research  does  not  have  an  actuator  modeled;  the  commanded  input  is 
directly  fed  through  to  the  model  without  being  augmented  by  actuator  dynamics. 
Therefore,  this  requirement  and  the  implications  associated  with  it  were  not  included  in 
the  scope  in  this  research. 

It  was  decided,  based  on  the  results  of  the  original  simulation  runs,  to  by-pass  the 
autopilot  in  the  simulation  with  a  command  generator  utilizing  a  series  of  different  pitch 
doublets  varying  in  frequency  including  band-limited  white  noise.  After  testing  multiple 
inputs  of  pitch  doublets  and  band-limited  white  noise  of  varying  frequencies  and 
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amplitudes,  two  input  types  showed  the  most  potential  to  be  used  for  the  estimation  of  the 
longitudinal  derivatives  of  the  missile.  These  input  types  were  chosen  because  they 
provided  the  most  excitation  to  the  system  without  exceeding  the  physical  limitations  of 
missile  actuator  hardware.  The  first  input  type  is  a  pitch  doublet  at  2  Hz,  shown  in  Figure 
4-3.  The  second  input  type  is  band-limited  white  noise  that  produces  a  normally 
distributed  random  signal,  shown  in  Figure  4-4. 


Pitch  Doublet  Input  Signal 


Figure  4-3  Two  Hz  Pitch  Doublet  Input  Signal 

White  Noise  Input  Signal 


8 
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Figure  4-4  Band-Limited  White  Noise  Input  Signal 
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4.3.2  Control  Signal  Input  Implementation  into  Simulation 

The  missile  is  modeled  in  the  cruciform  shape  (PHIPRM  =  0.0)  having  4  different 
control  surfaces  along  the  X  and  Y-axis  as  seen  in  Figure  4-5. 


PHIPRM  =  0.0  deg 


Figure  4-5  Cruciform  Missile  Layout 

The  relationship  between  the  individual  control  deflections  and  three  axes  control 
deflections  is  shown  from  Equation  3.10,  repeated  here. 

<51 1  ri  sin<zi  -cos^l,- 

Y  Y  SP 

52  1  cosd>  sin<z b 

V  V  SQ 

5  3  1  -sin<z b  cos<z J 

F  r  SR 

<54j  [1  -cos  (j)  -sin^J 

Fins  numbered  2  and  4  are  the  pitch  control  surfaces  and  Fins  1  and  3  control  the  roll 
and  yaw  axes.  Fins  2  and  4  are  coupled  together  with  Fin  4  commanding  opposite  that  of 
Fin  2.  This  equates  to  a  negative  command  input  for  Fin  4  in  the  simulation.  Since  this 
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research  focused  only  on  the  determining  the  longitudinal  stability  and  control 
derivatives,  Fins  number  2  and  4  were  commanded  with  the  pitch  command  inputs 
described  in  Figures  4-3  and  4-4  above  and  Fins  numbered  1  and  3  were  set  to  zero. 

4.4  Linear  Regression  Estimation 

The  data  was  collected  by  running  the  simulation  twice;  each  instance  using  one  of 
the  two  input  signals  selected  from  Section  4.3.  The  High  frequency  pitch  doublet  in 
Figure  4-3  will  be  labeled  Test  Case  1  and  the  Band  limited  white  noise  signal  in  Figure 
4-4  will  be  labeled  Test  Case  2.  The  variables,  a  (angle  of  attack  in  degrees),  q  (pitch 


rate  in  degrees/sec),  &  (Elevator  angle  in  degrees)  and  q  (pitch  acceleration  in 
degrees/sec")  were  collected  from  the  simulation  time  history  for  the  research.  They  are 
as  follows  for  Test  Case  1  and  Test  Case  2,  respectively. 
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Figure  4-6  High  Frequency  Pitch  Doublet  Input/Output  (Test  Case  1) 
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Figure  4-7  Band  Limited  White  Noise  Input/Output  (Test  Case  2) 

Data  collected  from  each  simulation  time  history  was  partitioned  into  discrete  samples 
taken  at  every  10th  sample,  or  every  0.01  seconds  for  the  entire  time  of  the  simulation, 
which  ran  for  10  seconds  with  10,000  samples  taken.  The  time  histories  were  partitioned 
in  this  way  for  the  estimation  of  the  longitudinal  dimensional  stability  derivatives,  Ma , 

Mq  and  . 

These  derivatives  are  computed  by  utilizing  the  Linear  regression  technique,  described  in 
Chapter  2,  to  solve  Equation  4. 1  from  Section  4.2 


q  =  Maa  +  Mqq  +  M &§e 

Using  Equation  2.9  from  Chapter  2  for  the  least-squares  estimator 

0  =  (XTXy1XTy 
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Where  X,  also  called  the  regressor  matrix,  and  Y  are  defined  as: 


Xl  =  \a  q  5e\ 


(4.13) 


Y  =  q  (4.14) 

Partitioning  the  data  to  every  10  samples  satisfies  the  requirement  that  m  >  n  in  order  to 
estimate  the  n  parameters  fort?  ,  where  m  =  10  and  n  =  3  .  From  Equation  2.4, 
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Figures  4-8  through  4-19,  which  will  be  discussed  in  the  following  section,  show  the 
comparison  between  the  estimated  values  of  Ma ,  Mq  and  M&  to  the  actual  values  from 

the  6-DoF  missile  model  simulation.  Some  of  the  Figures  show  the  effect  of  passing  the 
data  through  a  second-order  filter,  given  by: 


Filter  = 


co„ 


s  +  2£a>ns  +  con 


(4.16) 


Where, 


co„  -2k  -2 


and 


2 


(4.17) 


These  filtered  estimates  are  shown  in  Figures  4-11,4-15  and  4-19  for  Test  Case  1  and 
Figures  4-13,  4-17  and  4-21  for  Test  Case  2. 
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4.4.1  Linear  Regression  Estimation  Analysis  of  Results 


This  section  presents  an  analysis  of  the  results  obtained  from  applying  the  Linear 
Regression  method  to  estimate  the  longitudinal  stability  and  control  derivatives  from  the 
simulation  truth  data. 

4. 4. 1.1  Ma,  Pitch  angular  acceleration  per  unit  angle  of  attack: 

Figures  4-10  and  4-11  show  the  estimates  of  Ma  for  Test  Case  1  compared  to  the 

truth  data  from  the  simulation.  Once  the  estimate  data  was  filtered,  it  was  easier  to  see 
the  trend  of  the  data  for  each  estimate.  The  overall  trend  is  good,  however  there  is  some 
bias  in  the  estimate  data.  Figures  4-12  and  4-13  show  the  estimates  of  Ma  for  Test  Case 

2  compared  to  the  same  truth  data  from  the  simulation.  The  filtered  data  shown  in  Figure 
4-13  show  that  for  the  time  history  captured,  the  estimate  follows  well  with  the  trend  with 
very  little  bias.  While  the  estimate  doesn’t  exactly  track  the  oscillations  of  the  truth  data, 
the  overall  mean  values  compare  well,  based  on  a  visual  inspection  of  the  figure. 

Overall,  results  obtained  for  Ma  were  better  with  the  band-limited  white  noise  command 
signal  input  than  the  pitch  doublet. 

4.4. 1.2  M  ,  Pitch  angular  acceleration  per  unit  pitch  rate: 

Figures  4-14  and  4-15  show  the  estimates  of  M  for  Test  Case  1  compared  to  the 

truth  data  from  the  simulation.  Again,  it  was  easier  to  see  the  trend  of  the  data  using  the 
filtered  estimate.  The  estimate  had  an  average  bias  of  about  0.2  compared  to  the  truth 
data  and  had  some  erratic  data  points  in  the  estimate  near  end  of  the  time  history.  These 
erratic  data  points  at  the  end  are  cause  by  an  ill-condition  X  matrix  for  the  estimates  at 
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this  time  section.  This  will  be  further  discussed  later  in  this  report.  Figures  4-16  and  4- 
17  show  the  estimates  of  M  for  Test  Case  2  compared  to  the  same  truth  data  from  the 

simulation.  The  filtered  data  shown  in  Figure  4-17  show  that  for  the  time  history 
captured,  the  had  an  average  bias  of  about  0.2,  identical  to  Test  Case  1.  Overall,  both 
Test  Cases  failed  to  estimate  M  with  any  real  accuracy.  This  was  expected,  however, 

because  M  is  generally  difficult  to  estimate  with  parameter  ID  techniques. 

4. 4.1.3  M Se,  Pitch  angular  acceleration  per  unit  elevator  (pitch  deflection)  input: 

Figures  4-18  and  4-19  show  the  estimates  of  M &  for  Test  Case  1  compared  to  the 

truth  data  from  the  simulation.  From  the  filtered  estimate  data,  it  was  easier  to  see  the 
trend  of  the  data  for  each  estimate.  The  overall  trend  is  good;  there  is  little  overall  bias  in 
the  data  with  the  exception  to  the  erratic  data  points  found  in  the  data  near  the  beginning 
and  at  the  end  of  the  time  history.  These  erratic  data  points  result  from  ill-conditioned 
regressor  matrices,  X  ,  which  will  soon  be  discussed. 

Figures  4-20  and  4-21  show  the  estimates  of  M&  for  Test  Case  2  compared  to  the 

same  truth  data  from  the  simulation.  The  filtered  data  shown  in  Figure  4-21  shows  that 
for  the  time  history  captured,  the  estimate  follows  well  with  the  trend  of  the  truth  data. 
There  is,  however,  a  small  horizontal  shift  in  the  data  compared  to  the  un-filtered 
estimate  data  in  Figure  4-20.  This  shift  or  bias  is  caused  by  the  second  order  filter 
implemented  to  smooth  the  original  estimate  data.  If  the  truth  data  is  also  filtered  with 
the  second  order  filter  the  shift  goes  away.  This  is  seen  in  Figure  4-22.  Overall,  results 
obtained  for  M &  were  better  with  the  band-limited  white  noise  command  signal  input 
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than  the  pitch  doublet.  This  next  section  will  discuss  what  an  ill-conditioned  regressor 
matrix  is  and  how  it  is  detennined. 

4.4. 1.4  Condition  number  of  Regressor  Matrix,  Xt 

The  condition  number  of  a  matrix  provides  an  indication  of  the  sensitivity  of  the 
solution  of  a  system  of  linear  equations  to  errors  in  the  data.  The  2-norm  condition  of  a 
rectangular  matrix  is  the  ratio  of  the  largest  and  smallest  singular  values.  For  rectangular 
matrices  with  full  column  rank: 

X  e  IRmxn ,  rank(X)  =  n  =»  k(X)  =  (4.18) 

a,mn(x) 

It  gives  an  indication  of  the  accuracy  of  the  results  from  matrix  inversion  and  the 
linear  equation  solution.  Values  near  1  indicate  a  well-conditioned  matrix.  Large 
condition  numbers  indicate  a  nearly  singular  matrix  or  an  ill-conditioned  matrix.  Figures 
4-8  and  4-9  give  the  condition  number  of  the  estimates,  6X  for  the  regressor  matrix,  Xi  for 
both  Test  Case  1  and  Test  Case  2.  The  highest  and  lowest  condition  numbers  for  both 
commanded  signal  input  types  for  Xi  are  shown  in  Table  4-2  below. 
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Table  4-2:  Comparison  of  Maximum  and  Minimum  Condition  Numbers  for  Xi 


Regressor  Matrix,  Xi 

Highest  Condition  Number 

Lowest  Condition  Number 

Test  Case  1 

3e6 

3e3 

Test  Case  2 

5e3 

10 

time,  sec 


Figure  4-8  Condition  Number  of  Xi  for  Test  Case  1 
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Figure  4-9  Condition  Number  of  Xi  for  Test  Case 

From  observation  of  Figure  4-9,  high  condition  numbers  occur  at  the  beginning  and 
the  end  of  the  time  history  data  set.  This  has  a  direct  correlation  to  the  shape  of  the 
estimates  of  Ma,  M  and  M&  for  Test  Case  1  in  Figures  4-11,  4-15  and  4-19,  where  the 

erratic  data  points  in  the  beginning  and  end  of  the  time  history  exist.  This  indicates  that 
the  pitch  doublet  signal  used  for  Test  Case  1  was  not  set  to  a  high  enough  frequency  to 
excite  all  the  modes  of  the  system,  which  therefore,  produces  the  singularities  in  the 
Regressor  Matrix,  X.  These  singularities  also  occur  when  the  commanded  pitch  signal 
crosses  zero  in  the  simulation. 

Figure  4-9  shows  an  average  condition  number  for  Test  Case  2  for  the  entire  time 
history.  This  trend  could  translate  into  the  biases  seen  in  Figures  4-13,  4-17  and  4-21  for 


Condition  of  the  Regressor  matrix  (X) 
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the  estimates  for  Test  Case  2.  From  Figures  4-8  and  4-9,  it  is  clear  that  the  band-limited 
white  noise  signal  does  a  better  job  of  exciting  the  modes  of  the  system  than  the  pitch 
doublet  due  to  the  lower  condition  numbers  (fewer  singularities)  shown  in  Figure  4-9. 
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Comparison  of  Malpha  truth  to  Malpha  Linear  Regression 


Figure  4-10  Test  Case  1  Comparison  of  6xto  “Truth  Data”  for  Ma 


Comparison  of  Malpha  truth  to  Malpha  Linear  Regression  (Filtered) 


Figure  4-11  Filtered  Test  Case  1  Comparison  of  6xto  “Truth  Data”  for  Ma 
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Malpha  Malpha 


Comparison  of  Malpha  truth  to  Malpha  Linear  Regression 


Figure  4-12  Test  Case  2  Comparison  of  8X  to  “Truth  Data”  for  Ma 


Comparison  of  Malpha  truth  to  Malpha  Linear  Regression  (Filtered) 


Figure  4-13  Filtered  Test  Case  2  Comparison  of  6xto  “Truth  Data”  for  Ma 
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Comparison  of  Mq  truth  to  Mq  Linear  Regression 


Figure  4-14  Test  Case  1  Comparison  of  Ox  to  “Truth  Data”  for  M q 


Comparison  of  Mq  truth  to  Mq  Linear  Regression  (Filtered) 


Figure  4-15  Filtered  Test  Case  1  Comparison  of  9xto  “Truth  Data”  for  Mq 
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Comparison  of  Mq  truth  to  Mq  Linear  Regression 


Figure  4-16  Test  Case  2  Comparison  of  Ox  to  “Truth  Data”  for  M q 


Comparison  of  Mq  truth  to  Mq  Linear  Regression  (Filtered) 


Figure  4-17  Filtered  Case  2  Comparison  of  6X  to  “Truth  Data”  for  Mq 
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Comparison  of  Mdelq  truth  to  Mdelq  Linear  Regression 


Figure  4-18  Test  Case  1  Comparison  of  6xto  “Truth  Data”  for  M& 


Comparison  of  Mdelq  truth  to  Mdelq  Linear  Regression  (Filtered) 


Figure  4-19  Filtered  Test  Case  1  Comparison  of  0xto  “Truth  Data”  for  M & 
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Comparison  of  Mdelq  truth  to  Mdelq  Linear  Regression 


Figure  4-20  Test  Case  2  Comparison  of  6X  to  “Truth  Data”  for  M& 


Figure  4-21  Filtered  Test  Case  2  Comparison  of  6xto  “Truth  Data”  for  M& 
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Comparison  of  Mdelq  truth  (Filtered)  to  Mdelq  Linear  Regression  (Filtered) 


Figure  4-22  Comparison  of  Filtered  Est.  to  Filtered  Truth  Data 


4.5  Linear  Regression  Estimation  Validation 

The  next  step  of  the  research  focused  on  examining  the  validity  of  the  calculated 
estimates  for  Ma ,  Mq  and  M& .  To  test  this,  the  time  histories  of  the  estimates,  0X 
were  fed  into  the  linear  model  of  the  system,  from  Equation  4.3: 


a 

R 


X 

V 

a 

+ 

'V 

Ma 

Mq_ 
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Where: 


Za~~'  Ma  (4.19) 

c 

Z,=--M,  (4.20) 

c 

Z5  =  —  -Ms  (4.21) 

e  c 

Where,  /  =  0.4696  m  (missile  length)  and  c  =  0. 1280  m  (missile  diameter). 

This  linear  model  was  then  simulated  in  SIMULINK  in  discrete  time,  using  the  same 
command  inputs  as  for  the  6-DoF  missile  model  simulation  (from  Figures  4-3  and  4-4) 
and  the  output  data  was  recorded.  Once  again,  the  variables  a  (angle  of  attack  in 

degrees),  q  (pitch  rate  in  degrees/sec),  &  (Elevator  angle  in  degrees)  and  q  (pitch 
acceleration  in  degrees/sec“)  were  collected  from  the  simulation  time  history  for  the 
research.  These  time  histories  for  both  Test  Case  1  and  Test  Case  2  are  shown  below  in 
Figures  4-23  and  4-24. 
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Figure  4-23  Pitch  Doublet  Input/Output  for  Linear  Model  (Test  Case  1) 
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Figure  4-24  White  Noise  Input/Output  for  Linear  Model  (Test  Case  2) 
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The  derivatives,  Ma ,  M  and  M&  were  then  estimated  for  a  second  time  from  the  time 

history  outputs  of  the  linear  model  in  Equation  4.3  with  the  original  estimates  acting  as 
the  truth  states  of  the  model. 

The  original  estimates  0X  were  compared  to  the  new  estimates  of  the  linear  model,  02 
to  check  the  estimation  technique.  The  new  estimates  02  were  filtered  using  Equations 
4-14  and  4-15  and  were  compared  to  the  original  estimates  made  from  the  6-DoF  missile 
model  simulation.  Figures  4-27  through  4-32  shows  the  comparisons  made  between  the 
original  estimate,  9X  to  the  linear  model  estimate,  02  for  both  Test  Case  1  and  Test  Case 
2. 

4.5.1  Linear  Regression  Estimation  Validation  Analysis  of  Results 

This  section  presents  an  analysis  of  the  results  obtained  from  applying  the  Linear 
Regression  method  to  estimate  the  longitudinal  stability  and  control  derivatives  from  the 
simulation  truth  data. 

4. 5. 1.1  Ma,  Pitch  angular  acceleration  per  unit  angle  of  attack: 

Figure  4-27  shows  the  estimates,  02  of  Ma  for  Test  Case  1  compared  to  the 
original  estimates,  0] .  The  overall  trend  did  not  match  well  with  0] ,  and  there  were  two 
transient  noise  periods,  one  at  the  beginning  of  the  estimates  and  one  at  the  end.  Figure 
4-28  show  the  estimates  of  Ma  for  Test  Case  2  compared  to  the  original  estimates,  0X . 
The  overall  trend  was  good.  The  original  estimate  is  noisier  than  the  new  estimate. 
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Overall,  results  obtained  for  M a  were  better  with  the  band-limited  white  noise  command 
signal  input  than  the  pitch  doublet. 

4. 5. 1.2  M  ,  Pitch  angular  acceleration  per  unit  pitch  rate: 

Figure  4-29  shows  the  estimates,  #2  of  M q  for  Test  Case  1  compared  to  the 
original  estimates,  Gx .  The  new  estimate  has  some  transient  noise  at  the  beginning  and 
the  end  of  the  time  history,  but  in  the  middle  compares  better.  These  erratic  data  points 
causing  the  transient  noise  in  the  estimates  are  cause  by  an  ill-condition  X  matrix  for  the 
estimates  at  this  time  section.  This  will  be  further  discussed  later  in  this  report.  Figure  4- 
30  shows  the  estimates,  of  M q  for  Test  Case  2  compared  to  the  original  estimates,  Gx . 

Other  than  the  short  transient  noise  occurring  at  the  beginning  of  the  data,  the  overall 
trend  of  the  new  estimate  follows  the  original  estimate  very  well.  Surprisingly,  the 
estimates  made  for  M  for  both  Test  Cases  were  much  better  than  were  the  original 

estimates  to  the  truth  data. 

4.5.13  M&,  Pitch  angular  acceleration  per  unit  elevator  (pitch  deflection)  input: 

Figure  4-31  shows  the  estimates,  #2  of  M &  for  Test  Case  1  compared  to  the 
original  estimates,  Gx .  Once  again,  there  are  the  transient  noise  caused  by  erratic  data  in 
the  beginning  and  the  end  of  the  time  history  of  estimates.  This  phenomenon  will  be 
discussed  in  Section  4. 5. 1.4. 

Figure  4-32  shows  the  estimates,  02  of  M &  for  Test  Case  2  compared  to  the 
original  estimates,  Gx .  While  there  is  some  bias  between  the  estimates,  this  is  one  of  the 
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best  matching  comparisons.  The  bias  is  caused  by  the  second  order  filter  used  to  smooth 


the  estimate  data,  see  Section  4. 4. 1.3.  The  trend  is  identical  between  both  estimates  as 
are  any  fluctuations  in  the  data.  Overall,  once  again,  the  results  obtained  lor  M &  were 

better  with  the  band-limited  white  noise  command  signal  input  than  the  pitch  doublet. 
This  next  section  will  present  the  condition  number  analysis  made  for  the  Linear 
Regression  Validation  experiment. 

4.5. 1.4  Condition  number  of  Regressor  Matrix,  X2 

The  definition  of  the  condition  number  of  the  matrix  was  originally  defined  in 
Section  4. 4. 1.4.  Figures  4-25  and  4-26  give  the  condition  number  of  the  estimates,  62  for 
the  regressor  matrix,  X2  for  both  Test  Case  1  and  Test  Case  2.  The  highest  and  lowest 
condition  numbers  for  both  commanded  signal  input  types  for  X2  are  shown  in  Table  4-3 
below. 


Table  4-3:  Comparison  of  Maximum  and  Minimum  Condition  Numbers  for  X2 


Regressor  Matrix,  X2 

Highest  Condition  Number 

Lowest  Condition  Number 

Test  Case  1 

7e9 

lie3 

Test  Case  2 

6e7 

40 

Condition  of  the  Regressor  matrix  (X) 
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10 
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time,  sec 


Figure  4-25  Condition  Number  of  X2  for  Test  Case  1 


Condition  of  the  Regressor  matrix  (X) 
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Figure  4-26  Condition  Number  of  X2  for  Test  Case  2 


From  observation  of  Figure  4-25,  high  condition  numbers  occur  at  the  beginning  and 
the  end  of  the  time  history  data  set  similar  to  Figure  4-26.  This  again  has  a  direct 
correlation  to  the  shape  of  the  estimates  of  Ma,  M  and  M&  for  Test  Case  1  in  Figures 

4-27,  4-29  and  4-31,  where  the  erratic  data  points  in  the  beginning  and  end  of  the  time 
history  exist.  This  is  another  indication  that  the  pitch  doublet  signal  used  for  Test  Case  1 
was  not  set  to  a  high  enough  frequency  to  excite  all  the  modes  of  the  system,  which 
therefore,  produces  the  singularities  in  the  Regressor  Matrix,  X. 

Figure  4-26  shows  a  large  concentration  of  high  condition  numbers  at  the  beginning 
of  the  time  history  data.  This  trend  directly  correlates  to  the  cause  of  the  large  transient 
noise  data  points  seen  in  the  beginning  of  Figures  4-28,  and  4-30  for  the  estimates  for 
Test  Case  2.  From  Figures  4-8,  4-9,  4-25  and  4-26,  it  is  clear  that  the  band-limited  white 
noise  signal  does  a  better  job  of  exciting  the  modes  of  the  system  than  the  pitch  doublet 
due  to  the  lower  condition  numbers  (fewer  singularities). 

Generally  the  results  of  02  (the  longitudinal  derivatives  estimated  from  the  linear  model 
of  the  0X  estimates),  proved  to  be  an  acceptable  validation  of  the  linear  regression 
technique. 
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Mq  Estimate  Comparison 


Figure  4-29  Test  Case  1  Comparison  of  dxtod2ior  M 

Mq  Estimate  Comparison 


Figure  4-30  Test  Case  2  Comparison  of  6X  to  d2  for  M 
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5.  Conclusions  and  Recommendations 


5.1  Conclusions 

Using  Linear  Regression  as  a  parameter  identification  method  to  determine  the 
longitudinal  dimensional  stability  derivatives  of  a  tactical  missile  was  presented  in  this 
document.  The  linear  regression  technique  was  applied  and  compared  to  actual 
simulation  data  to  determine  the  suitability  of  the  method.  Two  control  input  forms  were 
tested  and  compared.  Estimates  were  determined  from  flight  test  data  obtained  from  a  6- 
DoF  missile  model  simulation  and  also  from  a  discrete  linear  model  of  the  system. 

This  research  shows  that  a  relatively  simple  estimation  method  such  as  linear 
regression  can  be  used  successfully  in  determining  longitudinal  dimensional  stability 
derivatives  of  a  tactical  missile  in  flight.  Comparative  results  presented  in  the  previous 
section  confirm  that  using  a  control  input  form  with  higher  frequency  modulations,  such 
as  band-limited  or  filtered  white  noise,  leads  to  good  estimates  of  model  parameters. 
While  this  method  is  not  as  complex  as  other  methods  of  parameter  estimation,  the  results 
show  that  satisfactory  estimates  can  be  obtained  in  a  short  amount  of  time.  This  could  be 
beneficial  for  engineers  who  require  a  rough  estimation  of  missile  aerodynamic 
parameters.  The  results  of  the  research  presented  here  will  help  shed  light  on  and  further 
the  study  of  parameter  identification  for  air-to-air  missiles. 

5.2  Recommendations  Further  Research 

The  linear  regression  method  presented  here  was  only  applied  to  determining  the 
longitudinal  dimensional  stability  derivatives  of  the  missile  simulation.  There  is  a  need, 
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however,  to  expand  this  research  to  include  the  lateral-directional  dimensional  stability 
derivatives  as  well. 

Future  work  should  also  include  improving  the  missile  simulation  used  in  this 
research.  The  missile  simulation  did  not  include  an  actuator  model  for  the  control 
surfaces  of  the  missile.  This  factor  could  significantly  alter  the  results  of  this  research. 
This  research  only  concentrated  on  data  input  control  signals  that  were  not  filtered  by 
actuator  dynamics.  The  actuator  model  might  introduce  physical  limitations  to  the  input 
signal,  which,  for  example,  could  limit  the  magnitude  and  frequency  of  the  input  signals 
that  were  used  in  this  research.  As  stated  earlier,  the  input  fonn  should  be  within  the 
bandwidth  of  the  actuator  driving  the  control  surface  being  commanded  for  the  maneuver 
of  the  missile.  This  factor  might  make  using  band-limited  white  noise  or  high  frequency 
pitch  doublets  as  input  signals  unfeasible  for  air-to-air  missile  parameter  identification. 

Expanding  the  experimental  factors  used  in  this  research  should  be  investigated. 
For  this  research,  only  two  forms  of  control  signal  actuation  were  examined.  More 
control  signal  input  types,  such  as  sinusoidal  inputs  or  chirp  signals  varying  in  frequency 
would  be  of  value  to  further  the  knowledge  on  this  subject  to  the  community.  Different 
sample  sizes  should  also  be  explored.  For  the  linear  regression  analysis,  samples  of  the 
time  history  data  sets  were  limited  to  only  10  samples  per  iteration,  or  every  0.01  seconds 
for  the  entire  time  of  the  simulation,  which  ran  for  10  seconds  with  10,000  samples  taken. 
The  best  estimates  resulting  from  this  research  were  for  the  derivatives,  Ma  and  M& . 

Further  research  needs  to  be  done  to  detennine  better  estimates  of  M  .  With  these 
recommendations  for  further  research  to  be  explored;  a  better  conclusion  can  be  made  on 
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the  applicability  of  using  Linear  Regression  analysis  for  tactical  missile  aerodynamic 
parameter  estimation. 
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